Elsevier

Computer Physics Communications

Volume 235, February 2019, Pages 169-178
Computer Physics Communications

Langevin dynamics simulation with dipole–dipole interactions: Massive performance improvements and advanced analytical integrator

https://doi.org/10.1016/j.cpc.2018.09.006Get rights and content

Abstract

The required computation duration is a crucial element for the computer simulation of both molecular and suspended nanoparticle dynamics. The latter one is already highly-optimized on a theoretical level by the implicit solvent notion within the Langevin and Brownian Dynamics techniques which still requires a more clear theoretical connection though. Also, the calculation performance issue becomes much more serious in case of the long-range interaction within the magnetostatics and/or electrostatics problems. Hence, both conceptual aspects and a technical challenge require a resolution for the real experimental applications. In order to do so, the analytical fluctuative–dissipative integrator has been derived from basic theoretical foundation. A numerical simulation complexity is partly reduced by the analytical solution. The method is reinforced by the Barnes–Hut general-purpose graphics processing unit (GPU) based dipole–dipole interaction algorithm. As a real application close to an experiment, the Langevin dynamics simulation of up to 1 million ferrofluid nanoparticles evolution have been performed at the all-with-all dipole–dipole long-range interactions without any cutoff radius and with real coefficients of translational and rotational viscous friction. Computational time scales comparisons of different methods have been made: 100 ns in silico could take between minutes till years of a computation Unix epoch time depending on a method. The 850 000 times calculation duration speed up from slowest (well-known) towards fastest (suggested here) simulation method has been achieved. Due to a GPU parallel computation the required calculation time scales up proportionally to a number of nanoparticles N which is much better compared to N*N or N*log(N) in case of a direct dipole summation or Barnes–Hut on the single processor method, respectively. The overall implementation is available and ready to be used within the wider open-source software package.

Introduction

Soft condensed matter and, particularly, complex fluid are of crucial interest and topicality in modern physical, chemical, and life sciences [[1], [2]]. They have also found applications in nanotechnology and medicine [[3], [4]]. Models of such matter are naturally based on statistical physics [[5], [6], [7], [8]]. Many practical calculations in the frameworks of such models are performed by means of a computer simulation: implementations of molecular dynamics (MD) [9], Monte Carlo [10], and hybrid Monte Carlo [11] algorithms. The partly modelled implicit solvent transforms MD concept into Langevin dynamics (LD) [[9], [12]], or to Brownian dynamics (BD) in case of the viscous limit approximation and getting rid of the inertial term of the stochastic differential equation of motion (Langevin equation) [13]. However, it is important to mention here that whereas MD and LD are numerical models of N,V,E-micro- and N,V,T-canonical ensembles, respectively, the BD does not fully model a statistical ensemble, rather its distinct aspect on the certain time scale. An introduction of a special friction term implying a mutual nanoparticles motion provides a Stokesian dynamics method [14]. A scaling need of the LD simulation towards a “mesoscopic” hydrodynamic description without the microscopic kinetic details loss had been implemented by dissipative particle dynamics (DPD) [[12], [15]] or other methods [16]. Any improvement of LD method can be easily generalized to DPD methods and other canonical ensemble models as it was done for a rotational motion quaternion based algorithm significant optimization [15]. The LD and DPD methods often are implemented within same simulation package and share many similar routines like Velocity Verlet (VV) integrator [[17], [18], [19]], interactions calculations, etc. Example of such package is “Extensible Simulation Package for the Research on Soft Matter Systems” (ESPResSo) [[20], [21], [22]].

Unlike BD, the LD simulation requires that the finite-difference scheme time-step order of magnitude should be smaller than the temporal parameters of the Langevin equation. Real physical problem (like Brownian motion of interacting nanoparticles) often corresponds to incommensurable characteristic times of a translational and rotational motions. As far as the time-step is common, this issue leads to an over-precision in one of these motion simulation and, finally, to a significant computational overhead. Probably, due to this technical problem, there are the research works which claim a relation between translational and rotational friction coefficients as a numerical method parameter for the equilibration time speed up. As a numerical approach of equilibria location, this approach is correct. Contrarily, a dynamic stochastic process consistent simulation, probably, should be performed based on real (i.e. derived from physical experiments or material parameters) translational and rotational friction coefficients, e.g. ones determined by Stokes’ law assuming the limit of low Reynolds number and the hydrodynamic-originated Langevin parameters [9]. Otherwise, numerical results could not be fully connected to real physical experiments: especially in case of dynamic processes investigations, electromagnetic radiation propagation and scattering, noise related phenomena, Brownian motion investigations, correct separation of metastable and the real equilibria, etc.

The method limitation of the incommensurable characteristic times requires a solution which will be suggested in the present work for ball dipole nanoparticles only. An LD simulation of the large ferrofluid aggregate is selected as a demo for such a challenge resolution. Such system as any other one with long-range each-with-each N nanoparticles’ interactions naturally corresponds to a calculation problem of a polynomial complexity of the ON2 type. The Barnes–Hut approximation [23] of the long-range dipole–dipole interaction implementation has been used in the present work.

Specific open-source MD package ESPResSo has been selected to publish a new method in a ready to be used way within a stable implementation with a large set of successfully passed automated validation tests scripted by other researchers from significantly different research programs within the computational soft matter science.

Section 2.1 is dedicated to the theoretical basis of a further derivation of known VV integrator and original one within Sections 2.2 Velocity Verlet fluctuative–dissipative integrator, 2.3 Analytical fluctuative–dissipative integrator , respectively. Section 3.1 describes parameters of real physical experiment applications: just to be as close to the physical reality as possible which is crucial in order to emphasize suggested analytical integrator capabilities for specific applied tasks. Section 3.2 reveals a technical implementation of the integrator in a way freely available for the open-source community. Finally, Section 4 is dedicated to a quantitative and illustrative demo of suggested integrator benefits and its critical comparison with other VV-like methods.

Section snippets

Theoretical foundation

The MD simulation method of the N,V,E-microcanonical ensemble should be symplectic (at an absence of fluctuative–dissipative components), energy-conservative [24], strictly time-reversible [18], and numerically stable [25]. These requirements are more simply achievable in practice for numerical integrators which are canonical transformations. Let us define the qj and the conjugate momentum pj as jth degrees of freedom canonical coordinates, whole set of which (system mechanical state) is

Materials and conditions

In order to test suggested AI, one is required to apply it to a real physical problem. One is also important to use material parameters and other experimental conditions which are similar to such reported by different researchers in the same applied science area. Let us consider a ferrofluid suspension [38] consisting of kerosene (dynamic viscosity η=1.64103 Pa s) as a carrier liquid and N magnetite (Fe3O4, a mass density ρ=5.24103 kgm3) ferromagnetic/superparamagnetic [39] nanoparticles

Results and discussion

The present simulation detailed setup parameters had been published in Tcl-script [53]. The original simplest cubic lattice structure of randomly switching dhS= 15.44 nm and dhL= 25.54 nm nanoparticles with the lattice parameter equal to the double large nanoparticle diameter have been selected in order to simulate the visual Stockmayer polymer-like structure quick enough formation [61]. The AI with/without current CUDA implementation of the BH method has been validated by the set of 111

Conclusions

Hence, a conjunction of Analytical fluctuative–dissipative integrator with GPU based Barnes–Hut long-range dipole–dipole interaction calculation simulation drastically increases Langevin simulation calculation capabilities: starting from 1 million nanoparticles a difference is minutes vs. years of Unix epoch time range. The Analytical fluctuative–dissipative integrator is a Velocity Verlet like modification of Brownian dynamics methods based on a formalism of the Liouville operator Trotter

Acknowledgements

We thank Mr. Alexander (Polyakov) Peletskyi and his scientific team for an approval of the Barnes–Hut GPU simulation source code publishing to open-source repositories. We thank Dr. Serhii Shulyma, Prof. Valeriy F. Kovalenko and Dr. Mykhaylo V. Petrychuk for experimental support, material parameters providing and discussions. We thank Mr. Nicola Tanygin for the simulation hardware appliance setup support. We thank ESPResSo and Ubuntu community for the open-source software maintenance which was

References (62)

  • FrenkelD.

    Physica A

    (2002)
  • JordanA. et al.

    J. Magn. Magn. Mater.

    (1999)
  • KolafaJ. et al.

    Fluid Phase Equilib.

    (1994)
  • LimbachH. et al.

    Comput. Phys. Comm.

    (2006)
  • TanyginB. et al.

    J. Magn. Magn. Mater.

    (2012)
  • HolmC. et al.

    Current Opinion in Colloid & Interface Science

    (2005)
  • ShulymaS. et al.

    J. Magn. Magn. Mater.

    (2016)
  • PolyakovA. et al.

    Comput. Phys. Comm.

    (2013)
  • IvanovA.O.

    J. Magn. Magn. Mater.

    (1996)
  • HamakerH.C.

    Physica

    (1937)
  • HumphreyW. et al.

    J. Mol. Graphics

    (1996)
  • JonesR.A.L.

    Soft Condensed Matter, Vol. 6

    (2002)
  • MashaghiS. et al.

    Int. J. Mol. Sci.

    (2013)
  • FisherI.Z.

    Statistical Theory of Liquids

    (1964)
  • FrenkelI.I.

    Kinetic Theory of Liquids

    (1955)
  • EconomouI.G.

    Ind. Eng. Chem. Res.

    (2001)
  • SchlickT.
  • PanagiotopoulosA.Z.

    Mol. Phys.

    (1987)
  • MehligB. et al.

    Phys. Rev. B

    (1992)
  • GogaN. et al.

    J. Chem. Theory Comput.

    (2012)
  • PottierN.

    Nonequilibrium Statistical Physics: Linear Irreversible Processes

    (2014)
  • BradyJ.F. et al.

    Annu. Rev. Fluid Mech.

    (1988)
  • MartysN.S. et al.

    Phys. Rev. E

    (1999)
  • O’ConnellS.T. et al.

    Phys. Rev. E

    (1995)
  • SwopeW.C. et al.

    J. Chem. Phys.

    (1982)
  • TuckermanM. et al.

    J. Chem. Phys.

    (1992)
  • OsaciM. et al.

    Microfluid. Nanofluid.

    (2017)
  • ArnoldA. et al.
  • The ESPResSo project Copyright (C), ESPResSo, 2016. URL:...
  • BarnesJ. et al.

    Nature

    (1986)
  • MurdockS.E. et al.

    J. Chem. Theory Comput.

    (2006)
  • Cited by (5)

    • Microscopic characteristics of magnetorheological fluids subjected to magnetic fields

      2020, Journal of Magnetism and Magnetic Materials
      Citation Excerpt :

      These steps iterate until the system reaches a stable equilibrium. In this paper, the velocity Verlet algorithm [36,37] is used to calculate the dynamical equation, which obtains the position and velocity of each particle simultaneously without losing accuracy. The velocity Verlet algorithm in this paper is stated as follows:

    • Microscopic characteristics of magnetorheological fluids for mining sampling devices in microgravity environment

      2021, Zhongnan Daxue Xuebao (Ziran Kexue Ban)/Journal of Central South University (Science and Technology)
    View full text