^†^†thanks: xuzl@sjtu.edu.cn^†^†thanks: zhaoyu14@msu.edu^†^†thanks: zhouqi1729@sjtu.edu.cn

Variance-reduced random batch Langevin dynamics

Zhenli Xu^1,2 Yue Zhao³ Qi Zhou¹ ¹School of Mathematical Sciences, Shanghai Jiao Tong University, Shanghai 200240, China
²CMA-Shanghai, MOE-LSC and Shanghai Center for Applied Mathematics, Shanghai Jiao Tong University, Shanghai 200240, China
³Department of Computational Mathematics, Science & Engineering, Michigan State University, East Lansing, MI 48824, USA

Abstract

The random batch method is advantageous in accelerating force calculations in particle simulations, but it poses a challenge of removing the artificial heating effect in application to the Langevin dynamics. We develop an approach to solve this issue by estimating the force variance, resulting in a variance-reduced random batch Langevin dynamics. Theoretical analysis shows the high-order local truncation error of the time step in the numerical discretization scheme, in consistent with the fluctuation-dissipation theorem. Numerical results indicate that the method can achieve a significant variance reduction since a smaller batch size provides accurate approximation, demonstrating the attractive feature of the variance-reduced random batch method for Langevin dynamics.

Langevin dynamics, random batch list method, variance reduction.

pacs:

02.70.Ns, 87.15.Aa

^†^†preprint: Preprint

I Introduction

Langevin dynamics simulation frenkel2001understanding is a highly influential computational technique utilized in various scientific disciplines, including physics, chemistry, biology, and material science, specifically at the nano- and micro-scales mabey2017strong ; Trullas1989Langevin ; Lee2017computational ; Pastor1994 . Its primary purpose is to investigate the equilibrium properties and dynamics of particles in many-body systems. As a thermostat in molecular dynamics simulation, this computational approach simulates canonical ensembles by integrating Newton’s equation of motion combined with thermal fluctuations and frictions welling2011bayesian ; ma2015complete ; roberts1996 ; dalalyan2017 ; oksendal03 , and its widespread applications have been found in the analysis of irreversible dynamics, stochastic properties, and trajectory instability of various systems norman2013stochastic . Furthermore, in cutting-edge fields such as machine learning, Langevin dynamics and its improved algorithms can efficiently transform stochastic optimization into posterior sampling, demonstrating strong potential for a wide range of learning tasks welling2011bayesian ; ma2015complete ; dalalyan2017 ; dubey2016variance ; coeurdoux2024normalizing .

The computational bottleneck for particle simulations using Langevin dynamics is non-bonded interactions between particles. The recently developed random batch method (RBM) jin2018random shows promise in the fast calculation of interacting particle systems. This method is inspired by the mini-batch idea in stochastic gradient descent algorithms of machine learning robbins1951stochastic ; bottou1998online ; bubeck2015convex , and approximates the force by using a batch of randomly selected neighbors of a fixed batch size such that the approximate force is unbiased, achieving much success CiCP-28-1907 ; li2020Random ; jin2022mean ; li2020some ; fornasier2022anisotropic ; ko2021model ; golse2019random ; ye2024error ; guillin2024some . Particularly, in molecular dynamics, the RBM has been extended to obtain the random batch Ewald method and the random batch sum-of-Gaussians method, which solve the parallel scalability for long-range interactions JinLiXuZhao2020 ; liang2022random ; liang2023random ; gan2024random . For short-range interactions, the singular kernel can be treated by the random batch list (RBL) method liang2021randomRBL ; liang2022IRBE ; lin2024hybrid . The RBL divides the interaction region into a core-shell structure and applies the RBM in the shell region, significantly reducing the neighboring number of particle interactions. The RBM is implemented for molecular dynamics at the graphic processing units and allows the simulation of a system of ten million particles in a single card gao2024rbmd .

The randomness introduced into the deterministic drifting term results in an additional error of $\mathcal{O}(\tau)$ at each discrete time step liang2021randomRBL , leading to the failure of the fluctuation-dissipation theorem kubo1966fluctuation . This can be well handled by introducing the Nosé-Hoover thermostat. However, under the Langevin thermostat, this stochastic approach can lead to an artificial heating effect. In the simulation of heterogeneous systems such as the gas-liquid coexistence, it often requires a larger batch size or decreasing time step size to maintain the equilibrium distribution of the system, reducing the efficiency of the RBM jin2022RBMSOI . The artificial heating effect raises the kinetic energy of particles to overcome the intermolecular forces, which may lead to unphysical outcomes in critical issues such as phase transition simulations guillin2024some .

In this paper, we present the variance-reduced random batch method for Langevin dynamics. The covariance matrix of the random forces in the drifting term is approximated and used to correct the Brownian motion, achieving a significant variance reduction and leading to the consistency of the single-step fluctuation-dissipation theorem. Theoretically, we demonstrate that the variance reduction method has higher discrete numerical accuracy in the statistical sense. The method is implemented into the Lennard-Jones fluid systems and Coulomb systems, and numerical results demonstrate the attractive performance with a significant reduction in the force variance. Notably, our algorithm exhibits linear computational complexity and effectively preserves the non-uniform properties observed in strongly correlated systems, thereby holding significant potential for large-scale simulations of heterogeneous systems.

II Method

Consider a system of $N$ particles at position $\bm{r}_{i}$ for $i=1,2,\cdots,N$ within a simulation box $\Omega$ with periodic boundary conditions. Suppose that the system is immersed in an infinite isothermal bath, and the particles obey the second-order Langevin dynamics,

m_{i}\mathrm{d}^{2}\bm{r}_{i}=\bm{F}_{i}\mathrm{d}t-\bm{\Gamma}\mathrm{d}\bm{r}_{i}+\bm{D}\mathrm{d}\bm{W}_{i}

(II.1)

where $\bm{r}_{i}$ is the particle position, $\bm{\Gamma}$ is the damping matrix, $\bm{D}$ is the diffusion matrix, and $\{\bm{W}_{i}\}$ is i.i.d. Brownian motion. The friction term $-\bm{\Gamma}\mathrm{d}\bm{r}_{i}$ and fluctuation force $\bm{D}\mathrm{d}\bm{W}_{i}$ satisfy the fluctuation-dissipation theorem kubo1966fluctuation ; marconi2008fluctuation ; prost2009generalized

2k_{B}T\bm{\Gamma}=\bm{D}^{2},

(II.2)

with $k_{B}T$ being the thermal energy. Introducing the particle velocity $\bm{v}_{i}$ with $\mathrm{d}\bm{r}_{i}=\bm{v}_{i}\mathrm{d}t$ , the classical Langevin dynamics Eq. (II.1) is discretized by the Euler-Maruyama scheme kloeden1992stochastic as follows,

\left\{\begin{array}[]{l}\Delta\bm{r}_{i}^{n}=\bm{v}_{i}^{n}\tau\\ m_{i}\Delta\bm{v}_{i}^{n}=\bm{F}_{i}^{n}\tau-\bm{\Gamma}\bm{v}_{i}^{n}\tau+\bm{\xi}_{i}^{n},\end{array}\right.

(II.3)

where $\Delta$ is the forward difference operator, $\tau$ is the time step, the superscript $n$ denotes the simulation step at $t=n\tau$ , and $\bm{\xi}_{i}^{n}\sim\mathcal{N}(0,\bm{D}^{2}\tau)$ denotes the discrete Brownian motion.

The total force acting on particle $i$ in Eq. (II.3) is composed of the external force and the interparticle forces,

\bm{F}_{i}=\bm{F}_{i}^{\mathrm{ext}}+\sum\bm{f}_{ij},

(II.4)

where one has $\bm{f}_{ij}=-\nabla_{\bm{r}_{i}}U_{ij}$ , and the summation runs over all particles in the central and image boxes. Typically, one considers the LJ potential

U_{ij}=4\epsilon_{\text{LJ}}\left[\left(\dfrac{\sigma_{\text{LJ}}}{r_{ij}}\right)^{12}-\left(\dfrac{\sigma_{\text{LJ}}}{r_{ij}}\right)^{6}\right],

(II.5)

where $\epsilon_{\text{LJ}}$ and $\sigma_{\text{LJ}}$ are the strength and the interaction distance scale. The classical Langevin dynamics introduces a cutoff radius $r_{s}$ , and particles beyond this distance are ignored in the force calculation. This results in a linear computational complexity by the linked cell list algorithm yao2004improved ; meloni2007efficient ; brown2011implementing ; welling2011efficiency , but for simulating the state of the gas-liquid coexistence, a large $r_{s}$ is required for accurate simulations. The RBL liang2021randomRBL has been proposed to reduce the number of interacting neighbors as well as the computational cost, by introducing the second cutoff radius $r_{c}$ , and dividing the interaction region into a core region $r<r_{c}$ , and a shell region $r_{c}\leq r<r_{s}$ . The force is directly calculated inside the core region, and the scaled force is computed in the shell zone by randomly selecting a batch of $P$ particles with uniform probability. Denote $N_{i}$ and $C_{i}$ as the particle number in the shell region and the set of the $P$ selected particles for the $i$ -th particle. The total force acting on this particle is given by

\widetilde{\bm{F}}_{i}=\bm{F}_{i}^{\mathrm{ext}}+\sum_{r_{ij}<r_{\mathrm{c}}}\bm{f}_{ij}+\frac{N_{i}}{P}\sum_{j\in C_{i}}\bm{f}_{ij},

(II.6)

and the average force is subtracted from each particle to conserve the total momentum.

In the RBL-based Langevin dynamics, the second equation of Eq. (II.3) reads

m_{i}\Delta\bm{v}_{i}^{n}=\widetilde{\bm{F}}_{i}^{n}\tau-\bm{\Gamma}\bm{v}_{i}^{n}\tau+\bm{\xi}_{i}^{n},

(II.7)

where the approximation Eq. (II.6) is used as the total force. It can be proved that the approximate force $\widetilde{\bm{F}}_{i}$ is unbiased with a bounded force variance and the convergent distribution is close to the stationary one liang2021randomRBL . The acceleration ratio can be approximated by the ratio of neighbor list size, as $[(4\pi/3)\rho r_{\mathrm{s}}^{3}]/[(4\pi/3)\rho r_{\mathrm{c}}^{3}+P]$ , where $\rho$ is the particle density.

Under the Euler-Maruyama scheme, the mean-field limit of the RBM results in an effective dynamics inass2023quantifying ; guillin2024some , where the random force introduces an additional contribution into the Brownian motion, such that $\bm{Z}_{t}\sim\mathcal{N}\left(0,\bm{D}^{2}\tau+\bm{\Sigma}\tau^{2}\right)$ . Here, $\bm{\Sigma}$ is the force covariance matrix. Therefore, it is necessary to reduce the force variance to maintain the fluctuation-dissipation theorem. One can approximately calculate the covariance matrix of the random force by

\begin{array}[]{c}\textbf{var}(\widetilde{\bm{F}}_{i})\approx\dfrac{N_{i}^{2}}{P^{2}}\left(\sum\limits_{j\in C_{P}}\bm{f}_{ij}\bm{f}_{ij}^{\mathrm{T}}-P\cdot\overline{\bm{F}}_{i}\overline{\bm{F}}_{i}^{\mathrm{T}}\right),\\ \overline{\bm{F}}_{i}:=\dfrac{1}{P}\sum\limits_{j\in C_{P}}\bm{f}_{ij}.\end{array}

(II.8)

Clearly, it is bounded due to the core-shell structure of the neighbor list and all particles in $C_{i}$ satisfying $r_{ij}\geq r_{c}$ . One subsequently subtracts the variance term from the Brownian motion in discrete Langevin dynamics with time step $\tau$ such that $\bm{\xi}_{i}^{n}$ in Eq. (II.7) is replaced by $\widetilde{\bm{\xi}}_{i}^{n}$ with

\widetilde{\bm{\xi}}_{i}^{n}\sim\mathcal{N}\left(0,\bm{D}^{2}\tau-\textbf{var}(\widetilde{\bm{F}}_{i}^{n})\tau^{2}\right),

(II.9)

where $\widetilde{\bm{F}}_{i}^{n}$ is calculated by (II.6) at the $n$ -th time step. The process of the Langevin dynamics evolved with the Brownian motion (II.9) is presented in Algorithm 1.

Algorithm 1 Varaince-reduced random batch Langevin dynamics

1: Choose

r_{s}

(the cutoff radius),

r_{c}

(the core radius),

\tau

(time step),

N_{t}

(total simulation steps), and

P

(batch size). Initialize positions, velocities and charges of all particles.

2: for

n\text{ in }1:N_{t}

3: Create the cell lists.

4: For each particle

i

, select

P

particles randomly into

C_{i}

from its shell zone if

N_{i}>P

, otherwise select all particles.

5: For each particle

i

, calculate the selected pair interactions

\bm{f}_{ij}

for all

r_{ij}<r_{c}

and

j\in C_{P}

j\neq i

, then calculate stochastic force

\widetilde{\bm{F}}_{i}

by (II.6) and covariance matrix

\textbf{var}(\widetilde{\bm{F}}_{i})

by (II.8).

6: Integrate the Langevin dynamics (II.7) with the Brownian motion (II.9).

7: end for

Considering the estimate of the local truncation error, one omits the superscript $n$ for simplicity. The velocity evolution of a single step in the RBL method is provided by (II.7), where $\widetilde{\bm{F}}_{i}$ is the unbiased estimation of $\bm{F}_{i}$ , hence it can be approximately regarded as $\widetilde{\bm{F}}_{i}=\bm{F}_{i}+\bm{\eta}_{i}$ with $\bm{\eta}_{i}\sim\mathcal{N}(0,\bm{\mathcal{V}})$ and $\bm{\mathcal{V}}$ denoting the covariance matrix of $\widetilde{\bm{F}}_{i}$ . The discrete scheme can be approximated as

m_{i}\Delta\bm{v}_{i}=(\bm{F}_{i}-\bm{\Gamma}\bm{v}_{i})\tau+\hat{\bm{\xi}_{i}},

(II.10)

where $\hat{\bm{\xi}_{i}}=\bm{\xi}_{i}+\bm{\eta}_{i}\tau+\mathcal{O}(\tau^{3/2})$ denotes the noise term. Hence the local truncation error of the RBL discretization is $\mathcal{O}(\tau)$ , which leads to the inconsistency with the fluctuation-dissipation theorem. As for the proposed dynamics with variance reduction technique, Theorem II.1 demonstrates that it reduces the local truncation error to $\mathcal{O}(\tau^{3/2})$ .

Theorem II.1.

Given the force variance, the Langevin dynamics evolved with the random force (II.8) and the Brownian motion (II.9) exhibits a local truncation error of $\mathcal{O}(\tau^{3/2})$ in the sense of statistical expectation.

Proof.

When applying the variance-reduced Brownian motion $\widetilde{\bm{\xi}}_{i}$ in (II.9) instead of $\bm{\xi}_{i}$ in the RBL dynamics (II.7), one finds that the noise term $\hat{\bm{\xi}}_{i}$ becomes

	$\displaystyle\hat{\bm{\xi}}_{i}$	$\displaystyle=\bm{\eta}_{i}\tau+\widetilde{\bm{\xi}}_{i}+\mathcal{O}(\tau^{3/2})$		(II.11)
		$\displaystyle=\bm{\gamma}_{i}+\mathcal{O}(\tau^{3/2}).$		(II.11)

where $\bm{\gamma}_{i}\sim\mathcal{N}(0,\bm{D}^{2}\tau+(\bm{\mathcal{V}}-\textbf{var}(\widetilde{\bm{F}}_{i}))\tau^{2})$ . The difference between $\bm{\gamma}_{i}$ and $\bm{\xi}_{i}$ is of $\mathcal{O}(\tau^{3/2})$ , since $\textbf{var}(\widetilde{\bm{F}}_{i})$ is an unbiased estimation of $\bm{\mathcal{V}}$ . So $\hat{\bm{\xi}}_{i}=\bm{\xi}_{i}+\mathcal{O}(\tau^{3/2})$ shows the higher order truncation error of the proposed dynamics. ∎

Along with insights from effective dynamics inass2023quantifying ; guillin2024some , one also obtains that the VR-RBL is in consistent with the fluctuation-dissipation theorem since it holds

2k_{B}T\bm{\Gamma}=\bm{D}^{2}+(\bm{\mathcal{V}}-\textbf{var}(\widetilde{\bm{F}}_{i}))\tau+\mathcal{O}(\tau^{3/2}).

(II.12)

Since the proposed dynamics have a higher order local truncation error, one can select a smaller batch size $P$ to achieve accurate simulation results. The primary computational complexity of variance-reduced dynamics centers around the force calculation of particles, while the computation of the variance reduction term entails comparatively lower computing overhead, requiring only the summation process of stored force arrays. The variance-reduced method is demonstrated to be accurate through numerical results of LJ fluid and colloid-ion systems in section III, indicating its potential as a general technique for the RBM that sheds new light on precise and effective Langevin dynamics for various applications.

III Results

We perform simulations to show the effectiveness of the proposed Langevin dynamics method by comparing the results of the three methods: classical Langevin dynamics (CL), the random batch list method without variance reduction (RBL), and the variance-reduced method (VR-RBL). The CL is a brute-force method and is supposed to be accurate as the reference solution. All the methods discretize the Langevin dynamics by the Euler-Maruyama scheme. We study two benchmark systems: one is the LJ fluid exhibiting a gas-liquid coexistence state, and the second is a colloid-ion system with varying correlation coefficients srinivas2004self ; beck2004methods ; miao2015accelerated . The simulations are conducted on a Linux system equipped with an Intel Xeon Scalable Cascade Lake processor running at 2.5 GHz, utilizing 1 single CPU core and 4 GB of memory.

Refer to caption — Figure 1: The RDF (a) and CPU time per step (b) calculated by the CL, RBL and VR-RBL methods. The reduced units are applied, and the lines in (b) indicate the linear fitting of all methods.

III.1 Lennard-Jones fluids

Consider the LJ fluid with interparticle potential given in (II.5) and parameters $\sigma_{\text{LJ}}=\epsilon_{\text{LJ}}=1$ . One sets the particle number $N=500$ in a cubic box such that the average particle density $\rho=0.2$ . The reduced units with $k_{B}=1$ , $\bm{\Gamma}=\bm{I}$ and time step $\tau=0.01$ are applied. The temperature takes $T=0.9$ at which the system exhibits gas-liquid coexistence state. The cutoff radius takes $r_{\mathrm{s}}=6.0$ such that the CL remains accurate at the gas-liquid coexistence state, and $r_{\mathrm{c}}=2.0$ . In the calculation, one performs $10^{5}$ time steps for equilibrium and then $2\times 10^{5}$ steps for statistics. Fig. 1(a) presents the radial distribution functions (RDF) by the three methods, where one calculates $P=5$ , $10$ , and $20$ for the RBL, and $P=5$ and $10$ for the VR-RBL. One can see that the variance reduction significantly improves the accuracy for given $P$ . The VR-RBL converges to the CL with the small batch size $P=5$ . In comparison, the RBL achieves the same accuracy at $P=20$ . In Fig. 1(b), we present CPU performance for the three methods by using the same setup with fixed density, but the number of particles increases up to $N=5\times 10^{6}$ . The computational cost of all methods increases linearly with the number of particles, but one can observe a significant reduction in the prefactor of the scaling law by imposing random batch techniques. One can roughly estimate that, for the LJ fluid system of $\rho=0.2$ , each particle has $125$ neighboring particles within $r_{\mathrm{s}}=6.0$ with $11$ neighbors in the core region. Fig. 1(b) shows that the RBL with $P=20$ is $4.6$ times faster than the CL, and for the VR-RBL the accelerating rate becomes $7.9$ . These results are in agreement with the rough estimate of neighboring number and demonstrate the high efficiency of the VR-RBL given the accuracy level.

III.2 Colloid-ion Systems

The second example is a colloid-ion system in a spherical simulation volume of radius $R=20$ . A colloidal particle with charge $Q_{0}=-20$ and radius $r_{0}=5$ is placed at the center, along with $100$ monovalent cations and $80$ monovalent anions with the same radius $r\equiv 0.5$ in the remaining volume. All quantities are provided in reduced units. The interparticle potential is composed of the LJ potential in (II.5) of parameters $\sigma_{\text{LJ}}=\epsilon_{\text{LJ}}=1$ , and Coulomb interaction $U_{ij}^{\mathrm{Coul}}=\alpha q_{i}q_{j}/r_{ij}$ with $\alpha$ being the coupling strength. A spherical Wigner-Seitz (WS) cell model groot1991ion ; tamashiro1998donnan is employed for the boundary condition in our simulation with the LJ for the ion-wall interaction. One sets the temperature $T=1.0$ , the unit damping matrix $\bm{\Gamma}=\bm{I}$ and time step $\tau=0.005$ in the CL, RBL and VR-RBL methods. In the simulation, one performs $5\times 10^{6}$ steps for equilibrium and another $5\times 10^{6}$ steps for statistics. The coupling strength takes $\alpha=2.0$ and $10.0$ , representing a weak and strong coupling system with corresponding cutoff radii setting $r_{\mathrm{c}}=4.0$ and $5.0$ , respectively. Due to the Coulomb interaction, the results of CL are derived from calculating all particle interactions.

Fig. 2 presents the results of the charge density of the cation (normalized by its bulk density) and the integrated charge distribution along the radial direction. The integrated charge distribution is calculated through $Q_{\mathrm{int}}(r)=Q_{0}+\sum\int_{0}^{r}q_{i}\rho_{i}(s)\mathrm{d}s$ , where $\rho_{i}(s)$ is calculated by a statistical average. One can observe that RBL and VR-RBL methods are accurate at the weak-coupling case of $\alpha=2.0$ (Fig. 2(a,c)). However, the RBL with a small batch size fails to capture the stationary distribution in the $\alpha=10.0$ case, which has strong Coulomb interactions (Fig. 2(b,d)). The VR-RBL with $P=5$ has almost the same accuracy as the RBL with $P=20$ , in agreement with the CL result. This result reveals the accuracy and efficiency advantages of the VR-RBL method.

It is also interesting to investigate how the time step affects the simulation accuracy. One considers the RBL and VR-RBL for the $\alpha=10.0$ case simulating the same times for the equilibrium and statistics, with the time step ranging from $\tau=0.001$ to $0.02$ . Fig. 3 presents the relative errors of the kinetic energies and potential energies in the system. These errors are calculated by statistical averages of the energies compared with those of the CL. One can clearly observe the linear convergence of the VR-RBL with $1/P$ for both the kinetic and potential energies. Both methods exhibit a convergence rate of $\mathcal{O}(\sqrt{\tau})$ ; however, the prefactor for the VR-RBL is significantly smaller than that of the RBL. Moreover, the VR-RBL method can accurately capture the property of the system with strong coupling strength even at larger $\tau$ , further highlighting its application potential.

IV Conclusions

In summary, we have proposed a variance-reduced random batch Langevin dynamics for the efficient and accurate simulations of particle systems. The VR-RBL computes interactions using the random batch method while simultaneously estimating the covariance matrix of the forces. The covariance matrix is used to correct the Brownian motion with little additional computational cost to mitigate artificial heating effects and maintain the consistency with the fluctuation-dissipation theorem. Numerical experiments conducted on LJ fluids with gas-liquid coexistence state and colloid-ion systems with different Coulomb strengths demonstrate the accuracy and computational efficiency of the VR-RBL. It is promising for simulations of large-scale heterogeneous systems, allowing for the use of larger time steps to reduce computational time. The rigorous analysis on the convergence of the VR-RBL, along with strategies for further optimizing its accuracy and computational efficiency, will be studied in our subsequent work. It is noteworthy that, the extension of our proposed dynamics to other random batch processes, such as the random batch Ewald and the random batch sum-of-Gaussians methods for long-range interaction systems in molecular dynamics simulations JinLiXuZhao2020 ; liang2023random , is straightforward to implement.

Acknowledgement

This work is supported by the National Natural Science Foundation of China (grant No. 12325113) and the Science and Technology Commission of Shanghai Municipality (grant Nos. 23JC1402300 and 21JC1403700). The authors also acknowledge the support from the HPC center of Shanghai Jiao Tong University.

Conflict of interest

The authors declare that they have no conflict of interest.

Data Availability Statement

The data that support the findings of this study are available from the corresponding author upon reasonable request.

References

(1) D. Frenkel, B. Smit, Understanding molecular simulation: From algorithms to applications, Vol. 1, Elsevier, 2001.
(2) P. Mabey, S. Richardson, T. White, L. Fletcher, S. Glenzer, N. Hartley, J. Vorberger, D. O. Gericke, G. Gregori, A strong diffusive ion mode in dense ionized matter predicted by langevin dynamics, Nature Communications 8 (1) (2017) 14125.
(3) J. Trullàs, A. Giró, J. A. Padró, Langevin dynamics simulation of ions in solution: Influence of the solvent structure, The Journal of Chemical Physics 91 (1) (1989) 539–545.
(4) J. Lee, H. Koh, D.-N. Kim, A computational model for structural dynamics and reconfiguration of dna assemblies, Nature Communications 14 (1) (2023) 7079.
(5) R. W. Pastor, Techniques and applications of Langevin dynamics simulations, Springer Netherlands, Dordrecht, 1994, pp. 85–138.
(6) M. Welling, Y. W. Teh, Bayesian learning via stochastic gradient Langevin dynamics, in: Proceedings of the 28th International Conference on Machine Learning (ICML-11), 2011, pp. 681–688.
(7) Y.-A. Ma, T. Chen, E. Fox, A complete recipe for stochastic gradient MCMC, in: Advances in Neural Information Processing Systems, Vol. 28, 2015, pp. 2917–2925.
(8) G. O. Roberts, R. L. Tweedie, Exponential convergence of Langevin distributions and their discrete approximations, Bernoulli 2 (4) (1996) 341–363.
(9) A. S. Dalalyan, Theoretical guarantees for approximate sampling from smooth and log-concave densities, Journal of the Royal Statistical Society Series B 79 (3) (2017) 651–676.
(10) B. Øksendal, Stochastic differential equations: An introduction with applications, 6th Edition, Springer, Berlin, Heidelberg, 2003.
(11) G. E. Norman, V. V. Stegailov, Stochastic theory of the classical molecular dynamics method, Mathematical Models and Computer Simulations 5 (2013) 305–333.
(12) K. A. Dubey, S. J Reddi, S. A. Williamson, B. Poczos, A. J. Smola, E. P. Xing, Variance reduction in stochastic gradient Langevin dynamics, Advances in Neural Information Processing Systems 29 (2016) 1154–1162.
(13) F. Coeurdoux, N. Dobigeon, P. Chainais, Normalizing flow sampling with langevin dynamics in the latent space, Machine Learning (2024) 1–26.
(14) S. Jin, L. Li, J.-G. Liu, Random batch methods (RBM) for interacting particle systems, Journal of Computational Physics 400 (1) (2020) 108877.
(15) H. Robbins, S. Monro, A stochastic approximation method, The Annals of Mathematical Statistics 22 (3) (1951) 400 – 407.
(16) L. Bottou, Online learning and stochastic approximations, in: Saad, D. (Ed.), On-line Learning in Neural Networks, Vol. 17, 1998, pp. 9–42.
(17) S. Bubeck, Convex optimization: Algorithms and complexity, Foundations and Trends in Machine Learning 8 (3-4) (2015) 231–358.
(18) S. Jin, X. Li, Random batch algorithms for quantum Monte Carlo simulations, Communications in Computational Physics 28 (5) (2020) 1907–1936.
(19) L. Li, Z. Xu, Y. Zhao, A random-batch Monte Carlo method for many-body systems with singular kernels, SIAM Journal on Scientific Computing 42 (3) (2020) A1486–A1509.
(20) S. Jin, L. Li, On the mean field limit of the random batch method for interacting particle systems, Science China Mathematics 65 (2022) 169–202.
(21) L. Li, J. G. Liu, Y. Tang, Some random batch particle methods for the Poisson-Nernst-Planck and Poisson-Boltzmann equations, Communications in Computational Physics 32 (1) (2022) 41–82.
(22) M. Fornasier, H. Huang, L. Pareschi, P. Sünnen, Anisotropic diffusion in consensus-based optimization on the sphere, SIAM Journal on Optimization 32 (3) (2022) 1984–2012.
(23) D. Ko, E. Zuazua, Model predictive control with random batch methods for a guiding problem, Mathematical Models and Methods in Applied Sciences 31 (08) (2021) 1569–1592.
(24) F. Golse, S. Jin, T. Paul, The random batch method for $n$ -body quantum dynamics, Journal of Computational Mathematics 39 (6) (2021) 897–922.
(25) X. Ye, Z. Zhou, Error analysis of time-discrete random batch method for interacting particle systems and associated mean-field limits, IMA Journal of Numerical Analysis 44 (3) (2024) 1660–1698.
(26) A. Guillin, P. Le Bris, P. Monmarché, Some remarks on the effect of the random batch method on phase transition, Stochastic Processes and their Applications 179 (2025) 104498.
(27) S. Jin, L. Li, Z. Xu, Y. Zhao, A random batch Ewald method for particle systems with Coulomb interactions, SIAM Journal on Scientific Computing 43 (4) (2021) B937–B960.
(28) J. Liang, P. Tan, L. Hong, S. Jin, Z. Xu, L. Li, A random batch Ewald method for charged particles in the isothermal–isobaric ensemble, The Journal of Chemical Physics 157 (14) (2022) 144102.
(29) J. Liang, Z. Xu, Q. Zhou, Random batch sum-of-Gaussians method for molecular dynamics simulations of particle systems, SIAM Journal on Scientific Computing 45 (5) (2023) B591–B617.
(30) Z. Gan, X. Gao, J. Liang, Z. Xu, Random batch Ewald method for dielectrically confined Coulomb systems, arXiv preprint arXiv:2405.06333.
(31) J. Liang, Z. Xu, Y. Zhao, Random-batch list algorithm for short-range molecular dynamics simulations, The Journal of Chemical Physics 155 (4) (2021) 044108.
(32) J. Liang, Z. Xu, Y. Zhao, Improved random batch Ewald method in molecular dynamics simulations, The Journal of Physical Chemistry A 126 (22) (2022) 3583–3593.
(33) H. Lin, Y. Shi, S. Dai, Hybrid random batch idea and nonlinear conjugate gradient method for accelerating charged polymer dynamics simulation, Journal of Mathematical Chemistry 62 (3) (2024) 555–578.
(34) W. Gao, T. Zhao, Y. Guo, J. Liang, H. Liu, M. Luo, Z. Luo, W. Qin, Y. Wang, Q. Zhou, S. Jin, Z. Xu, RBMD: A molecular dynamics package enabling to simulate 10 million all-atom particles in a single graphics processing unit, arXiv preprint arXiv:2407.09315.
(35) R. Kubo, The fluctuation-dissipation theorem, Reports on Progress in Physics 29 (1) (1966) 255.
(36) S. Jin, L. Li, Y. Sun, On the random batch method for second order interacting particle systems, Multiscale Modeling & Simulation 20 (2) (2022) 741–768.
(37) U. M. B. Marconi, A. Puglisi, L. Rondoni, A. Vulpiani, Fluctuation–dissipation: Response theory in statistical physics, Physics Reports 461 (4-6) (2008) 111–195.
(38) J. Prost, J.-F. Joanny, J. M. Parrondo, Generalized fluctuation-dissipation theorem for steady-state systems, Physics Review Letters 103 (9) (2009) 090601.
(39) P. E. Kloeden, E. Platen, Stochastic differential equations, in: Numerical Solution of Stochastic Differential Equations, Springer, 1992, pp. 103–160.
(40) Z. Yao, J.-S. Wang, G.-R. Liu, M. Cheng, Improved neighbor list algorithm in molecular simulations using cell decomposition and data sorting method, Computer Physics Communications 161 (1-2) (2004) 27–35.
(41) S. Meloni, M. Rosati, L. Colombo, Efficient particle labeling in atomistic simulations, The Journal of Chemical Physics 126 (12) (2007) 121102.
(42) W. M. Brown, P. Wang, S. J. Plimpton, A. N. Tharrington, Implementing molecular dynamics on hybrid high performance computers – short range forces, Computer Physics Communications 182 (4) (2011) 898–911.
(43) U. Welling, G. Germano, Efficiency of linked cell algorithms, Computer Physics Communications 182 (3) (2011) 611–615.
(44) S. Inass, G. Stoltz, Quantifying the mini-batching error in Bayesian inference for adaptive Langevin dynamics, The Journal of Machine Learning Research 24 (1) (2024) 15638–15695.
(45) G. Srinivas, D. E. Discher, M. L. Klein, Self-assembly and properties of diblock copolymers by coarse-grain molecular dynamics, Nature Materials 3 (9) (2004) 638–644.
(46) D. A. Beck, V. Daggett, Methods for molecular dynamics simulations of protein folding/unfolding in solution, Methods 34 (1) (2004) 112–120.
(47) Y. Miao, F. Feixas, C. Eun, J. A. McCammon, Accelerated molecular dynamics simulations of protein folding, Journal of Computational Chemistry 36 (20) (2015) 1536–1549.
(48) R. D. Groot, Ion condensation on solid particles: Theory and simulations, The Journal of Chemical Physics 95 (12) (1991) 9191–9203.
(49) M. N. Tamashiro, Y. Levin, M. C. Barbosa, Donnan equilibrium and the osmotic pressure of charged colloidal lattices, The European Physical Journal B-Condensed Matter and Complex Systems 1 (1998) 337–343.