The Basics of Molecular Simulations: Part-1
The fundamental idea behind molecular simulations is that every system in the universe is composed of atoms connected to each other. These atoms interact with others which can be described with certain equations and with the help of these equations, the system properties can be obtained. Molecular simulation can be applied for a wide variety of systems such as biological, materials, etc. Below are some common applications,
- In the field of material science, material properties such as density, volume, conductivity, etc. change based on external conditions. Molecular simulations can predict these changes.
- In the field of drug discovery, molecular simulations are used to study antibody-antigen reactions to identify drug candidates.
Generally, molecular simulations are used to complement experimental findings and gain molecular-level insights that are not possible experimentally.
Molecular simulations must be carried out in such a manner that they resemble the external conditions in which real phenomena/experiment is happening. For that matter, ensembles are defined. The most common ensembles are Isothermal-isobaric (N, P and T kept constant), canonical (N, V, and T kept constant), microcanonical (N, V, and E kept constant), and grand canonical (μ, V and T kept constant).
N = Total number of molecules in the system
P = Total Pressure of the system
T = Temperature of the system
V = Volume of the system
E = Total energy of the system
μ = Chemical potential of the system
Before explaining further, I would like to explain my research in short and how I used molecular simulations for my research. I was researching on adsorption of CO2 in Metal-organic frameworks (MOFs). MOFs are porous materials, found to have excellent adsorption properties. I used Grand Canonical Monte Carlo (GCMC) and Molecular Dynamics (MD) simulations to understand the system consisting of MOF and CO2 atoms.
Grand Canonical Monte Carlo (GCMC): GCMC simulation technologies are useful when we want to obtain equilibrium properties of the system of interest for e.g. adsorption capacity, vapor-liquid phase equilibrium properties, etc. In GCMC simulation, we specify the temperature and chemical potential of the species of interest, for a fixed simulation volume (μ, V and T kept constant). GCMC works on the Markov chain method. Markov model is a stochastic model describing a sequence of possible events in which the probability of each event depends only on the state attained in the previous event. The event is defined as the introduction of the new atoms inside the system. The probability of moving to the next configuration is proportional to the energies of the respective configurations.
In figure 1.1, two configuration N and N+1 are shown. N+1 configuration was obtained after inserting one atom at a random position in configuration N. The probability of moving to configuration N+1 from N, is calculated based on the energies of these two configurations, and successively the system will evolve [1].
I used GCMC to calculate sorption isotherms. These isotherms provide us the amount of adsorbate (X) adsorbed on the adsorbent surface with different pressures (P) at a constant temperature (T). At each pressure, the equilibrium configuration is obtained by adding some molecules to the configuration at the previous pressure. It gives us the adsorption capacity at that particular pressure as depicted in figure 1.2.
Molecular Dynamics (MD): MD simulations are used when we want to understand the dynamic behavior of the system. Atoms’ motion will be predicted following the classical newton’s law of motion and successively the change in the system properties takes place. MD methods are suitable for those in which we want to study the dynamic properties of the system such as diffusivity, conductivity, etc. They can also be used to compute the equilibrium properties but GCMC is much faster as it doesn’t require solving complex equations and just works on calculating the probability. The potential energy between atoms is calculates using the formula written below,
Where;
Ubond = oscillations about the equilibrium bond length
Uangle = oscillations of 3 atoms about an equilibrium bond angle
Udihedral = torsional rotation of 4 atoms about a central bond
Unonbond = non-bonded energy terms (electrostatics and Lenard-Jones)
The interactions between atoms can be classified into two categories, bonded and non-bonded. Van der Waals and electrostatic attractions are non-bonded interactions. Bonded interactions are the movements of atoms around bonds. For a system of N particles with coordinates X and velocities V, the following pair of first-order differential equations may be written in Newton’s notation as
The potential energy function U(X) of the system is a function of the particle coordinates X. It is referred to simply as the force field. The first equation comes from Newton’s laws of motion; the force F acting on each particle in the system can be calculated as the negative gradient of U(X). For every time step, each particle’s position X and velocity V are updated. The time evolution of X and V is called a trajectory. Given the initial positions (e.g., from theoretical knowledge) and velocities (e.g., randomized Gaussian), we can calculate all future (or past) positions and velocities. From these values, we can acquire system properties. The thermal conductivity calculation of argon employing MD simulations in the NVT ensemble is the most studied example [2].
Summary:
Figure 1.3 below summarizes both MC and MD simulations techniques;
References:
In the next part, I will cover the tools which I used to perform molecular simulations and talk more about the parameters needed for them.
