Metropolis Monte Carlo (MMC) Algorithm
The MMC algorithm implements sequential steps that allow the chosen atomic sphere to move or not to a new spatial position consistent with the atomic bonding structure. The list of the (x, y, z) positions for the centers of each individual sphere were updated after each MMC step as these are the values to be used for the desired calculations of average thickness of a monolayer and the formation of gauche-bonds. Using classical physics, the algorithm computes the total energy of a sphere at a particular “initial” position and also at a randomly chosen new position, a “final” position. The algorithm specifies the calculations to be done to either accept or reject the possible new-final position.
The total energy of a sphere is given by: i) the bond stretching energy,\(E_{b}\), between neighboring bonded spheres (on the same molecule) and, ii) the dispersion energy, \(E_{d}\), between non-bonded spheres on the same molecule, or spheres on different molecules . Thus, the total energy \(E,\ \) for a particular sphere is
\(E=\ E_{b}+\ E_{d}\) [1]
Where \(E_{b}\), was represented by
\(E_{b}={3.8\ \ \left(r-3\right)}^{4}\text{\ \ }\) [2]
Here the unit is one-half Å, and r represents the center-to-center distance of bonded (adjacent) spheres. The energy of 3.8 eV is the energy of a single C-C bond, . The value of 3.8 (eV/units4) was used to cause the bonding energy to reach 3.8 eV at a distance of ± 1 unit from the minimum.
\(E_{d}\) was calculated from equation 3
\(E_{d}=\left(\frac{j}{r^{8}}-\frac{k}{r^{6}}+l\right)\frac{5}{l}\ \)[3]
Where, as above, r represents the distance between two non-bonded spheres. r6 is the well-known attractive London dispersion energy . The short-range repulsion ofr-8 was used rather than the typicalr-12 (. We chose r-8 in order to permit the formation of gauche bonds, by making the electron cloud less repulsive than when using r-12 ─ if no significant number of gauche bonds arise usingr-8 , then they will not arise usingr-12 . \(j,\ k\ \text{and}\text{\ l}\) are coefficients to assign an energy of 0 as the minimum, based on the equilibrium distance of 4.9 Å between parallel hydrocarbon chains .
Notice that the total energy does not include electrostatic energy, as this was ignored so that our simulation excludes BL.
In this system, the solvent oil is assumed to be a neutral background which does not play a significant role in the formation of gauche bonds, hence it was not necessary to consider it.
An MMC algorithm starts by selecting a sphere to move. This simulation followed a sequential order starting from the “center” sphere in each molecule. This center of the molecule was held constant throughout the simulation. Molecules were prevented from moving perpendicular to the long axis (Fig 1).
PLACE FIGURE 1 HERE
Since each molecule’s center was held fixed throughout the simulation, its location and the distances between atoms, and therefore, molecules were known.
An MMC step was started by recording the spatial position of the selected sphere (initial position) and by computing its initial energy\(\ E_{i}\) (eq. 1). Note that the \(E_{b}\ \)(eq. 2) component was computed using the spheres bonded to the selected sphere in this initial position, while \(E_{d}\) (eq. 3) was computed using the “non-bonding” neighboring spheres. A random number was then used to select a possible new position (final position) for the selected sphere. The final energy \(E_{f}\) was computed for this group at the new coordinates. The newly computed \(E_{f}\) was subtracted from \(E_{i}\)(eq. 4)
\(\Delta E=E_{f}-E_{i}\) [4]
The final calculation in the MMC step was to determine if the new position was to be accepted by looking at the energy change \(\Delta E\)(eq. 4) value. If \(\Delta E\) was negative, a lower energy state was achieved, hence the move to the new position was accepted. If the change in energy was positive, a new calculation was made to decide if this higher energy state was acceptable. Equation 5, which takes into account the entropy of the system, details the calculation made
\(g=\ e^{(\frac{-\Delta E}{k_{B}T})}\ \) [5]
where kB is Boltzmann’s constant (eV) and T = 300 °K for this simulation.
A random number R was then chosen in the range \(0\leq R<1\). If \(R\leq g\), then the new position was accepted if the move was allowed. However, if R > g the new position was not accepted.
After the last calculation, a new sphere was selected and all the previous steps repeated. Spheres were selected sequentially starting at the center. Once all spheres in the molecule were treated with one MMC step, the software routine comes back to the first central-fixed group and started over. Each group was attempted to be moved sufficient times for a total of 100,000 times (105 MMC steps). Ten replicates were run and averaged to obtain the results.