i-PI: A Python interface for ab initio path integral molecular dynamics simulations

doi:10.1016/j.cpc.2013.10.027

Computer Physics Communications

Volume 185, Issue 3, March 2014, Pages 1019-1026

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

Abstract

Recent developments in path integral methodology have significantly reduced the computational expense of including quantum mechanical effects in the nuclear motion in ab initio molecular dynamics simulations. However, the implementation of these developments requires a considerable programming effort, which has hindered their adoption. Here we describe i-PI, an interface written in Python that has been designed to minimise the effort required to bring state-of-the-art path integral techniques to an electronic structure program. While it is best suited to first principles calculations and path integral molecular dynamics, i-PI can also be used to perform classical molecular dynamics simulations, and can just as easily be interfaced with an empirical forcefield code. To give just one example of the many potential applications of the interface, we use it in conjunction with the CP2K electronic structure package to showcase the importance of nuclear quantum effects in high-pressure water.

Program summary

Program title: i-PI

Catalogue identifier: AERN_v1_0

Program summary URL: http://cpc.cs.qub.ac.uk/summaries/AERN_v1_0.html

Program obtainable from: CPC Program Library, Queen’s University, Belfast, N. Ireland

Licensing provisions: GNU General Public License, version 3

No. of lines in distributed program, including test data, etc.: 138626

No. of bytes in distributed program, including test data, etc.: 3128618

Distribution format: tar.gz

Programming language: Python.

Computer: Multiple architectures.

Operating system: Linux, Mac OSX, Windows.

RAM: Less than 256 Mb

Classification: 7.7.

External routines: NumPy

Nature of problem:

Bringing the latest developments in the modelling of nuclear quantum effects with path integral molecular dynamics to ab initio electronic structure programs with minimal implementational effort.

Solution method:

State-of-the-art path integral molecular dynamics techniques are implemented in a Python interface. Any electronic structure code can be patched to receive the atomic coordinates from the Python interface, and to return the forces and energy that are used to integrate the equations of motion.

Restrictions:

This code only deals with distinguishable particles. It does not include fermonic or bosonic exchanges between equivalent nuclei, which can become important at very low temperatures.

Running time:

Depends dramatically on the nature of the simulation being performed. A few minutes for short tests with empirical force fields, up to several weeks for production calculations with ab initio forces. The examples provided with the code run in less than an hour.

Introduction

Molecular dynamics simulations are becoming increasingly capable not only of assisting the interpretation of experiments, but also of predicting the properties and the behaviour of new materials and compounds [1]. Besides the increase in available computer power, these developments have been made possible by an increasingly accurate treatment of the interactions between the atoms, in particular by an explicit treatment of the electronic structure problem [2], [3]. However, as the methods that are used to treat the quantum nature of the electrons improve, it becomes increasingly clear that in the presence of light atoms, such as hydrogen or lithium, the error due to approximating the nuclei as classical particles is at least as large as the errors due to the approximate modelling of the electrons. The importance of nuclear quantum effects (NQEs) is evident from the fact that the zero-point energy associated with a typical O–H stretching mode is in excess of 200 meV. This has significant implications: for example, one can show by extrapolating the experimental values for isotopically pure light, heavy and tritiated water that a hypothetical liquid with classical nuclei would have a 50% higher heat capacity than H₂O and a pH of about 8.5.

Within the Born–Oppenheimer approximation, NQEs can be modelled using the imaginary time path integral formalism [4], [5], [6], [7]. This formalism maps the quantum mechanical partition function for a set of distinguishable nuclei onto the classical partition function of a so-called ring polymer, composed of $n$ replicas (beads) of the physical system connected by harmonic springs. The number of replicas needed to achieve a converged result is typically a small multiple of $β ħ ω_{max}$ , where $β$ is the inverse temperature and $ω_{max}$ is the largest vibrational frequency in the system. For the archetypical case of room-temperature liquid water, most properties are reasonably well converged with $n = 32$ —making the path integral simulation 32 times more expensive than a classical simulation.

As a consequence of this large overhead, NQEs were for many years only rarely considered in the context of ab initio molecular dynamics [8], [9], [10]. However, with the advent of massively parallel computers, including these effects has become somewhat more affordable [11], [12], and their simulation has also been facilitated by new methodological developments. In particular, it has recently been realised that an approximate modelling of NQEs can be obtained by applying a coloured (correlated)-noise Langevin equation thermostat to classical molecular dynamics [13], [14], [15]. Moreover, this idea can be turned into a systematically convergent, accurate method [16], [17] by combining coloured noise with path integral molecular dynamics (PIMD). To give an example, the PIGLET method [17] makes it possible to treat NQEs in liquid water at 300 K with as few as 6 beads, which promises to make accurate ab initio PIMD calculations almost routine.

Unfortunately, most ab initio electronic structure codes only contain rudimentary implementations of PIMD, if any at all. We have developed i-PI in order to reduce the effort of introducing the latest PIMD developments into electronic structure codes, so as to make them readily available to a broader community. The framework we have adopted is a server–client model, in which i-PI acts as a server passing atomic coordinates to (multiple instances of) an electronic structure client program, and receiving the energy and forces in return. This allows all of the PIMD machinery to be confined to the Python server side, and minimises the modifications that have to be made to the electronic structure code.

Section snippets

Program overview

i-PI has been developed with the awareness that the cost of an ab initio PIMD simulation will be dominated by the electronic structure problem, and with a few clear goals:

•
to minimise the effort required to modify the client ab initio electronic structure code,
•
to include no feature specific to a particular client code,
•
to be modular, simple to extend and with a structure that reflects the underlying physics,
•
to be efficient, but never at the cost of clarity.

In order to achieve these goals we have

Program features

Integrating the PIMD equations of motion involves a number of technicalities. Fortunately, these can be kept to a minimum by using stochastic thermostatting [18] — which avoids the complication of integrating Nosé–Hoover chain thermostats — and by formulating the time evolution using a symmetric Trotter splitting algorithm [19]. This is the approach we have adopted in the i-PI program.

The program implements the most recent developments in path integral and coloured-noise molecular dynamics,

Constant-pressure path integral molecular dynamics

Most of the techniques listed above have been discussed in detail elsewhere. However, our implementation of PIMD in the NPT ensemble has not been described before now. The approach we have adopted for this is a relatively straightforward combination of ideas taken from path integral Langevin equation thermostats [18], stochastic barostats for conventional MD [34], and previous work on constant-pressure PIMD [35]. We report it here as this combination is robust, transparent and streamlined,

An example application: high-pressure water

To highlight some of the more advanced features available in i-PI, we have used it to perform an illustrative ab initio PIGLET [17] simulation of supercritical water at 750 K and 10 GPa. This simulation is both technically challenging and physically interesting. We shall use it to reveal the importance of NQEs under conditions similar to those explored in the pioneering work of Ref. [39], in which calculations were performed at constant volume and without quantum effects in the nuclear motion.

Conclusions

In this paper we have introduced i-PI, a Python interface designed to facilitate including nuclear quantum effects in ab initio path integral molecular dynamics simulations. Our program delegates the calculation of the potential, forces and virial tensor to an external code, keeping the electronic structure calculation and the propagation of the nuclear dynamics separate. Communication between the codes is achieved using internet sockets, which exchange just the essential information, thereby

Acknowledgements

We would like to thank the early adopters of i-PI, including J. Cuny, E. Davidson, R. DiStasio, D. Donadio, F. Giberti, A. Hassanali, T. Markland, M. Rossi, B. Santra, D. Selassie, T. Spura, L. Wang, and D. Wilkins, who have helped immensely debugging the code and testing it on a variety of problems. A special thanks goes to A. Michaelides, who helped us to decide an acronym for the interface. We also acknowledge generous allocations of computer time from CSCS (project ID s338) and the Oxford

References (52)

T.E. Markland et al.
A refined ring polymer contraction scheme for systems with electrostatic interactions
Chem. Phys. Lett.
(2008)
M. Ceriotti et al.
The $δ$ -thermostat: selective normal-modes excitation by colored-noise Langevin dynamics
Procedia Comput. Sci.
(2010)
C.C. Fischer et al.
Predicting crystal structure by merging data mining with quantum mechanics
Nat. Mater.
(2006)
R. Car et al.
Unified approach for molecular dynamics and density-functional theory
Phys. Rev. Lett.
(1985)
K. Burke
Perspective on density functional theory
J. Chem. Phys.
(2012)
D. Chandler et al.
Exploiting the isomorphism between quantum theory and classical statistical mechanics of polyatomic fluids
J. Chem. Phys.
(1981)
M. Parrinello et al.
Study of an F center in molten KCl
J. Chem. Phys.
(1984)
D.M. Ceperley
Path integrals in the theory of condensed helium
Rev. Modern Phys.
(1995)
R.P. Feynman et al.
Quantum Mechanics and Path Integrals
(1964)
M.E. Tuckerman et al.
Efficient and general algorithms for path integral Car{–}Parrinello molecular dynamics
J. Chem. Phys.
(1996)

GLE4MD

Cited by (253)

Revisiting first principles van der Waals corrections based on maximally localized Wannier functions
2024, Chemical Physics Letters
The dependence of structural and dynamical properties of ambient water on a newly introduced scaling factor $γ_{s}$ for the formulation of van der Waals (vdW) corrections based on maximally localized Wannier functions was systematically investigated via conventional first principles MD and ring-polymer MD at constant NPT conditions. We find that this small modification improves the description of both the properties. This occurs along with a slight reduction of the density to $\sim$ 0.97 g/cm³, provided that the value of $γ_{s}$ is appropriately selected. This correction allows for a significant improvement toward a reliable and quantitative performance of first principles simulations.
Hydrogen trapping and embrittlement in metals – A review
2024, International Journal of Hydrogen Energy
Hydrogen embrittlement in metals (HE) is a serious challenge for the use of high strength materials in engineering practice and a major barrier to the use of hydrogen for global decarbonization. Here we describe the factors and variables that determine HE susceptibility and provide an overview of the latest understanding of HE mechanisms. We discuss hydrogen uptake and how it can be managed. We summarize hydrogen trapping and the techniques used for its characterization. We also review literature that argues that hydrogen trapping can be used to decrease HE susceptibility. We discuss the future research that is required to advance the understanding of HE and hydrogen trapping and to develop HE-resistant alloys.
CRYSTALpytools: A Python infrastructure for the CRYSTAL code
2023, Computer Physics Communications
CRYSTALpytools is an open source Python project available on GitHub that implements a user-friendly interface to the Crystal code for quantum-mechanical condensed matter simulations. CRYSTALpytools provides functionalities to: i) write and read Crystal input and output files for a range of calculations (single-point, electronic structure, geometry optimization, harmonic and quasi-harmonic lattice dynamics, elastic tensor evaluation, topological analysis of the electron density, electron transport, and others); ii) extract relevant information; iii) create workflows; iv) post-process computed quantities, and v) plot results in a variety of styles for rapid and precise visual analysis. Furthermore, CRYSTALpytools allows the user to translate Crystal objects (the central data structure of the project) to and from the Structure and Atoms objects of the pymatgen and ASE libraries, respectively. These tools can be used to create, manipulate and visualise complicated structures and write them efficiently to Crystal input files. Jupyter Notebooks have also been developed for the less Python savvy users to guide them in the use of CRYSTALpytools through a user-friendly graphical interface with predefined workflows to complete different specific tasks.
Program Title: CRYSTALpytools
CPC Library link to program files: https://doi.org/10.17632/p2bp3fsk86.1
Developer's repository link: https://github.com/crystal-code-tools/CRYSTALpytools
Licensing provisions: MIT
Programming language: Python and Jupyter Notebook
Nature of problem: The Crystal code [1,2] is a powerful tool for the calculation of materials properties. It stands out in the condensed matter computational landscape because of the use of local basis sets, heavy parallelisation, efficient implementation of non-local Fock exchange, and extensive use of point and space symmetry. However, it currently lacks an easily programmable interface to access its input/output structure needed to be able to use Crystal calculations within computational materials science workflows. Historically, the use of Crystal in such workflows has been achieved through bespoke scripting implemented in a variety of languages which has hindered code reuse and sharing.
Solution method: The CRYSTALpytools project will enable the automation of Crystal calculations in a modular code that can be co-developed by a wide community of users worldwide. It achieves this by transforming standardised input and output files into python objects and providing a suite of functionality to manipulate them. The core implementation is based on a set of data structures denoted the Crystal_objects. CRYSTALpytools, in its current implementation, contains a large variety of functions for input/output manipulation, vibrational and thermodynamic analysis, and visualisation.
Additional comments including restrictions and unusual features: The exploitation of point group and space group symmetry is one of the strengths of the Crystal code, therefore, all the geometry modification functions have been developed to optimise the use of symmetry. When a structure is downloaded from a database or is modified by external code, such as pymatgen [3] or ASE [4], a symmetry analysis is performed before transforming the structure to be an optimised Crystal_object. All such objects can be directly transformed to standard structures periodic in 0, 1, 2 or 3 dimensions used by Crystal (ie: .gui or .f34 files).
- [1]
  Dovesi Roberto, et al., Quantum-mechanical condensed matter simulations with CRYSTAL, WIREs Comput. Mol. Sci. 8(4) (1 July 2018) e1360.
- [2]
  Erba Alessandro, et al., CRYSTAL23: A Program for Computational Solid State Physics and Chemistry, J. Chem. Theory Comput., 2022.
- [3]
  Shyue Ping Ong, et al., Python Materials Genomics (pymatgen): A Robust, Open-Source Python Library for Materials Analysis, Computational Materials Science 68 (2013) 314-319.
- [4]
  Ask Hjorth Larsen, et al., The atomic simulation environment—a Python library for working with atoms, J. Phys.: Condens. Matter, 29 (7 June 2017) 273002.
Path-integral molecular dynamics predictions of equilibrium H and O isotope fractionations between brucite and water
2023, Geochimica et Cosmochimica Acta
Precise prediction of hydrogen (H) isotope fractionations among different substances is a long-standing challenge in isotope geochemistry, as it needs treatments beyond the harmonic approximation. Path-integral molecular dynamics (PIMD) simulations have recently been proved to be valid in predicting equilibrium isotope fractionations for light elements. However, the lack of reliable force fields hinders the application of PIMD to the condensed phases. In this work, the deep potential models trained on the first-principles molecular dynamics (FPMD) data are applied to PIMD simulations to determine the D/H and ¹⁸O/¹⁶O isotope fractionation factors between brucite and water. To quantitatively assess the influence of quantum effects, the D/H and ¹⁸O/¹⁶O isotope fractionations between brucite and water are also determined under the framework of harmonic approximation.
By comparing the results of the PIMD and those of harmonic calculations, we find that H and O atoms in brucite and water are strongly affected by the anharmonic effects. The accuracy of the harmonic isotope fractionations mainly depends on the cancellation percentage of the anharmonic effects on the reduced partition function ratios (RPFRs) of the two substances. The isotope fractionations predicted by PIMD simulations are close to the experimental results at high temperatures. However, at low temperatures, the predicted isotope fractionations are different from those of experiments. The discrepancies are attributed to the approximations in the density functional theory (DFT) functionals used in this work, i.e., the inaccuracy in predicting the proton positions in brucite and water. Moreover, we find that the H isotope fractionations are extremely sensitive to the change of pressure. The direction of H isotope fractionation between brucite and water (at 300 K) will be even inversed as the pressure increases to more than 3 GPa. This strong pressure sensitivity may be a common characteristic of hydrous minerals and water systems. Therefore, the extrapolated H isotope fractionations used in the studies of dehydration of subducting slabs or deep Earth H distribution based on the isotope fractionations under low pressures may need to be rechecked.
Fortnet, a software package for training Behler-Parrinello neural networks
2023, Computer Physics Communications
A new, open source, parallel, stand-alone software package (Fortnet) has been developed, which implements Behler-Parrinello neural networks. It covers the entire workflow from feature generation to the evaluation of generated potentials, coupled with higher-level analysis such as the analytic calculation of atomic forces. The functionality of the software package is demonstrated by driving the training for the fitted correction functions of the density functional tight binding (DFTB) method, which are commonly used to compensate the inaccuracies resulting from the DFTB approximations to the Kohn-Sham Hamiltonian. The usual two-body form of those correction functions limits the transferability of the parametrizations between very different structural environments. The recently introduced DFTB+ANN approach strives to lift these limitations by combining DFTB with a near-sighted artificial neural network (ANN). After investigating various approaches, we have found the combination of DFTB with an ANN acting on-top of some baseline correction functions (delta learning) the most accurate one. It allowed to introduce many-body corrections on top of two-body parametrizations, while excellent transferability to chemical environments with deviating energetics could be demonstrated.
Program title: Fortnet
CPC Library link to program files: https://doi.org/10.17632/sjg3n9vr8p.1
Developer's repository link: https://github.com/vanderhe/fortnet
Code Ocean capsule: https://codeocean.com/capsule/3992747
Licensing provisions: LGPL
Programming language: Fortran, Python
External routines/libraries: MPI, BLAS/LAPACK, HDF5, DFTB+
Supplementary material: See supplementary material for exemplary Human-friendly Structured Data (HSD) input listings, as well as the basic usage of the Fortformat Python layer for generating datasets and extracting results.
Nature of problem: Semi-empirical quantum mechanical methods like density functional tight binding (DFTB) rely on fitting empirical energy correction terms, often represented by two-body potentials, to ab initio references. Hereby empirical, beyond-pairwise contributions are inevitably incorporated and therefore inadequately covered by a purely two-body description.
Solution method: The new, open source, parallel, stand-alone software package Fortnet provides a powerful, yet accessible tool to construct many-body correction terms by resorting to high-dimensional neural networks of Behler-Parrinello type. Fortnet is characterized by its modern infrastructure, complementing the landscape of available implementations by a robust combination of Fortran and Python based code.
Additional comments including restrictions and unusual features: Fortnet's core is supplemented by two additional projects that are BSD 2-clause licensed, namely fortnet-python [2], a collection of Python based tools for generating compatible datasets and extracting results, and fortnet-ase [3], an interface to the Atomic Simulation Environment (ASE) [1]. Both projects are available via the Python Package Index (PyPI). The interaction of all components is explained in cookbook-like recipes (see: https://fortnet.readthedocs.io/en/latest/), meant to guide new users, while learning about various basic and more advanced features by using comprehensible examples with physical reference.
- [1]
  A.H. Larsen, J.J. Mortensen, J. Blomqvist, I.E. Castelli, R. Christensen, M. Dułak, J. Friis, M.N. Groves, B. Hammer, C. Hargus, E.D. Hermes, P.C. Jennings, P.B. Jensen, J. Kermode, J.R. Kitchin, E.L. Kolsbjerg, J. Kubal, K. Kaasbjerg, S. Lysgaard, J.B. Maronsson, T. Maxson, T. Olsen, L. Pastewka, A. Peterson, C. Rostgaard, J. Schiøtz, O. Schütt, M. Strange, K.S. Thygesen, T. Vegge, L. Vilhelmsen, M. Walter, Z. Zeng, K.W. Jacobsen, J. Phys., Condens. Matter 29 (2017) 273002. doi:10.1088/1361-648x/aa680e.
- [2]
  T. van der Heide, fortnet-python, (2022), GitHub repository, https://github.com/vanderhe/fortnet-python.
- [3]
  T. van der Heide, fortnet-ase, (2022), GitHub repository, https://github.com/vanderhe/fortnet-ase.
Three-Center Tight-Binding Together with Multipolar Auxiliary Functions
2024, Journal of Chemical Theory and Computation

View all citing articles on Scopus

^☆: This paper and its associated computer program are available via the Computer Physics Communication homepage on ScienceDirect (http://www.sciencedirect.com/science/journal/00104655).

View full text

i-PI: A Python interface for ab initio path integral molecular dynamics simulations☆

Abstract

Program summary

Introduction

Section snippets

Program overview

Program features

Constant-pressure path integral molecular dynamics

An example application: high-pressure water

Conclusions

Acknowledgements

Chem. Phys. Lett.

Procedia Comput. Sci.

Predicting crystal structure by merging data mining with quantum mechanics

Nat. Mater.

Unified approach for molecular dynamics and density-functional theory

Phys. Rev. Lett.

Perspective on density functional theory

J. Chem. Phys.

Exploiting the isomorphism between quantum theory and classical statistical mechanics of polyatomic fluids

J. Chem. Phys.

Study of an F center in molten KCl

J. Chem. Phys.

Path integrals in the theory of condensed helium

Rev. Modern Phys.

Quantum Mechanics and Path Integrals

Efficient and general algorithms for path integral Car{–}Parrinello molecular dynamics

J. Chem. Phys.

Ab initio path integral molecular dynamics: Basic ideas

J. Chem. Phys.

The nature of the hydrated excess proton in water

Nature

Proton transport in triflic acid hydrates studied via path integral car-parrinello molecular dynamics

J. Phys. Chem. B

Nature of proton transport in a water-filled carbon nanotube and in liquid water

PCCP

Computation of the temperature dependence of the heat capacity of complex molecular systems using random color noise

Phys. Rev. E

Nuclear quantum effects in solids using a colored-noise thermostat

Phys. Rev. Lett.

Quantum thermal bath for molecular dynamics simulation

Phys. Rev. Lett.

Accelerating the convergence of path integral dynamics with a generalized Langevin equation

J. Chem. Phys.

Efficient first-principles calculation of the quantum kinetic energy and momentum distribution of nuclei

Phys. Rev. Lett.

Efficient stochastic thermostatting of path integral molecular dynamics

J. Chem. Phys.

Reversible multiple time scale molecular dynamics

J. Chem. Phys.

An efficient ring polymer contraction scheme for imaginary time path integral simulations

J. Chem. Phys.

Canonical sampling through velocity rescaling

J. Chem. Phys.

Langevin equation with colored noise for constant-temperature molecular dynamics simulations

Phys. Rev. Lett.

Colored-noise thermostats à la Carte

J. Chem. Theory Comput.