Low Hanging Fruit in Computational Molecular Dynamics Simulations with the Published OPLS-AA Force Field

After performing computations in molecular dynamics in high performance research systems, many solvable but yet unsolved problems have been determined. These include (i) a list of molecules which are of commercial and scholarly interest that can be parameterized with the available OPLS-AA literature, (ii) a list of metals in certain environments which have yet to be parameterized and would provide insight into their behavior in commercial processes, (iii) an empirical approach to searching for the impact of molecular parameterizations on observable properties which is useful in searching for new materials, and (iv) understanding molecular reactivity in computational molecular dynamics. *Corresponding author: Mongelli GF, Department of Chemical Engineering, Case Western Reserve University, 10900 Euclid Ave., A.W. Smith 116, Cleveland, USA, Tel: 12163682000; E-mail: Gfm12@case.edu Received December 14, 2016; Accepted December 30, 2016; Published January 10, 2017 Citation: Mongelli GF (2016) Low Hanging Fruit in Computational Molecular Dynamics Simulations with the Published OPLS-AA Force Field. J Material Sci Eng 6: 312. doi: 10.4172/2169-0022.1000312 Copyright: © 2017 Mongelli GF. This is an open-access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.


Introduction
Recently, the use of computational molecular dynamics packages, GROMACS in particular, has opened many avenues for innovative studies into chemical properties of materials [1]. Several studies by Jorgensen [2], Sides [3] and Frischnecht and Curro [4] have pointed out the key parameters or molecular potentials. These potentials are employed together with theory and high performance computing resources to produce ensemble average molecular properties which, often times mirror experimentally derived ones [5]. Combining the software with the potentials creates a basis set of permutations of these parameterizations, which represent all the simulatable materials.
Researchers have yet to probe the ensemble averages of a small fraction, let alone the entirety of the available set. Neither have scholars listed or prioritized which potentials are key to theoretical materials studies. This manuscript addresses both of those considerations. It highlights several key, classical materials which are within the span of such a basis parameterization set.
After spending some time within a particular subfield, typically a researcher may be able to list many of the solvable problems they can see combining the known input parameterizations, the software tools available and the computational power available.
Considering polymeric materials in the context of the OPLS-AA force field gives insight into their surfactant capability. The guideline for solubility in pure aqueous systems is that if the vector sum of the dipole moments or a molecular structure has an ensemble average that is zero, then such a molecule would exhibit hydrophobicity in water. The addition of an alcohol co-solvent solvating polyether materials in most cases polymers is a result of the 'like dissolving like' formalism. The carbon components of the alcohols have favorable non-bonded interactions with the polymers such as poly-alkanes, poly-ethers and poly-silicone which facilitate their solubilization.
Many problems associated with molecular simulations in organic or polymer systems are associated with system scale. Some of these problems that would be beneficial to solve remain unapproached due to their conceptual complexity and numerical scale. Over time, these will become cheaper to solve in accordance with Moore's law, however, each presents its own unique scientific challenges. Correspondingly, each provides a more complete understanding of how statistical mechanics broadly describes macromolecular behavior.

Discussion
Demonstrating the capability to correctly evaluate the surface activity opens up a new avenue for us to perform arbitrary structure searches for new materials with desirable properties. It remains to be shown if the ensemble average surface tensions reproduce experiments for systems large enough to make variations in the pressure tensor negligible.
There are some key deviations our theoretical work has had from industrial systems. For example, the author is aware of no simulations of surfactants in multi-solvent systems include sanitizing materials which may alter the surface activity and surface tensions in a concentrationdependent way. One such sanitizing material of interest is isopropyl myristate. This material is simulateable within the OPLS-AA force field. A vitamin E additive commonly added to hand sanitizers, tocopheryl acetate, is of interest. Tocopheryl acetate is simulateable within the OPLS-AA force field. These molecules can and should be incorporated into molecular dynamics simulations.
Cases where a polymer which are soluble in organic solvent A and solvent B but not their mixture AB can be searched for in future studies. The author is not aware of any such case. Solvophobic polymers in pure solvents A and B which become solvophilic in solvent AB are of interest as well. These could have unique applications in the surfactants industry and is related to packing and space filling in hydrogen bonding networks and would necessitate discontinuities in the surface excesses or free energy of solvation as a function of mass content. Given what is known about discontinuities in specific heats of particular liquid crystal mixtures, and broken ergodicity, the possibilities of such phenomena may not be impossibly surprising.
The inputs for the polyether and polysilicone materials are included and can be readily simulated in these solvent mixtures. The surface parameter of these polymers at higher molecular weights could also be studied. It is possible to scale to polymers to large molecular weight through a simple recursive generation program, yet to be written. Another tool exists, rdb files with pdb2gmx, however, it is poorly documented.

Force field determination tasks and gaps
GROMACS software makes an approximation in utilization of LJ potentials and lacks the complete multi-pole description of van der Walls forces. It should be considered which ensemble properties are most heavily determined by which terms such that more accurate computations can be performed when desired. More than this, those properties can be linked to their atomistic forces. This may have implications in the prediction of the exact form of liquid crystal phases and modes.
Additionally, there are many missing transition metal force field potentials within the OPLS-AA series. Several metals of note include Ir, Mo, Ag, Au, Pt, and Pd. These metals are present in many organic, biologically and optically active materials. Furthermore, C=Si systems could be parameterized. Determination of such potentials will offer new insights into the behavior of light emitting materials and inorganicorganic complexes. The parameterization of such metal complexes and alloys associated forces will allow for the determination of which systems will take best advantage of the Cashmir effect and relativistic retardations in extraction of work from the zero point energies of systems [6].
There are many related materials of interest to simulate properties of such as silanes, especially those with analogies to PE, PP, and PIB. It might yield interesting insights to study the surface tension reduction capabilities of materials and solvent mixtures with varied steric strength. Examples of such side groups are methyls, alcohols, carbonyls and styrene groups. Potential solvents include un-simulated hydroxyacetone and other double bond systems in addition to dimethylether. The dimethyl ether can be considered the poly(ethyleneoxide) with the lowest possible degree of polymerization.
A complex and rigorous method to determine the limitations of solubilization of molecules is to systematically turn off the potentials within each of the soluble and insoluble molecular structures and determine the ensemble average properties of such molecules. Such a method tells us immediately the potential surfactantability of compounds associated atoms in newly synthesized environments as soon as their OPLS-AA force field constants are determined experimentally. An alternative is to scale their force constants spatial dependence through the constant or not varying the strength at all but simply scaling back a certain type of force on an atom type or molecule. Then, to break the differences in force field parameters into intervals between those ranges for C and O and Si, and shift one molecular structures topology onto another in individual and all potentials. In effect, determining the dose response of molecular interactions on their macromolecular properties. This will proffer the contributions of various interactions to the discernable properties and give great insight into molecular design. Unfortunately, the computational resources and code structures to perform such an analysis have not been developed and are several years away at best.
Studies which change the dihedral constants while maintaining the degree of polymerization allow for the study of the fraction of the configurational space, and therefore number of configurations, achievable by a certain number of backbone atoms relative to the solubility. In such a manner, the entropy contributions to the solubility effect from the dihedral interactions can be separated from the Lennard Jones and Coulombics.
Alternatively, to alter the Lennard Jones constants while holding the Coulombic constants unchanged and studying the same tendencies. Such a study would essentially break down the mixing of spatial energy drop-offs of unfavorable and favorable interactions contributing to the ensemble average properties. It would be of principal interest to see how radial distribution function, free energy of solvation and surface activity vary when such force field parameters are altered to give rise to theoretical molecules.

Needed simulation automation tools
Through improvements to the OPLS-AA force field and better software links between computational molecular dynamics methods and quantum force field determination derivations, the inputs for molecular dynamics simulations could be readily achievable for every element in the periodic table and every molecule in an organic or inorganic chemistry text that is a liquid at STP. Most notable and interesting to the author are polymeric materials, drug compounds, particularly those whose binding sites are unknown within the body, and organic light emitting materials.
Future studies could look at the surface activity of known odor causing chemicals to determine if use of potentially newly determined surfactant structures would result in mixtures that can still maintain their cleansing properties and solubilization of antibacterial molecules. Some of the odor causing chemicals as determined by the Preti group [7] is butyric acid, trans-3-methyl-2-hexanoic avid, 5-hexenoic acid, androstenone, and pentyl acetate. Some potential solvents of interest include hydroxyacetone, methylglyoxal, peruvaldehude, methylglyoxal 1,1,-dimethyl acetal and glyoxal.

Modified surface active quantitative tools
Future FORTRAN codes to write include: the fraction of the interfacial slab which contains each of the two solvents as derived from density profiles or directly from coordinates of atoms. Also, to determine if there are any spatial dependencies of the hydrogen bonding -i.e., within 1.5 nm of interfaces --to see if more hydrogen bonding occurs at the interface or in the bulk. Writing a software program to compare the surface activity of a molecule across alcohol contents and pair that to the concentration of each solvent as a function of z would give insight into spatial molecular partitioning and how that scales with surface area of a particular solvent mixture and its size. This alludes to the concept of surface excess free energy.

Insight for materials synthesis heuristics
The trends in aqueous solubility and solvophilicity analogies for polyalkanes and polyethers could be tested. In other words that certain numbers of carbons or oxygens between silicon predict surface activity and surface tension reduction. Furthermore that certain side groups in such structures trigger surface tension reduction. If they fail to do so, new insight may be gained regarding molecular design. Within the polysilicone series, some of the required dihedral interaction class parings associated with analogues to common polyether and polyalkane surfactants do not have known energy potential parameterizations. Quantum molecular dynamics packages may be employed to determine such potentials.
Un-substituted saturated alkane polymers are ideal for comparison with polyethers. Such polymers are topologically similar to and can be mapped onto can be readily heteroatom substituted -or mapped -to become topologically similar to polyethers and polysilicones. The choice of the alkane structures is such as to study the varied branching, molecular weight, and methyl side-group density of carbon structure and molecular weight.

More automated scalable and increasingly complex systems
Potential future studies include the evaluation of multiple polymer molecule and polydisperse systems. Polydispersity would be interesting because of the opportunity to determine its impact on simulated viscoelastic macromolecular behavior, in particular dynamic mechanical analyses. Other simulations could vary the acidity/basicity and ionic character of solutions or seek to determine the minimum number of polymer molecules at a certain chain length to form a stable micelle. The surface activity of these polymeric materials in mixed solvents such as butanols, propylene glycol, glycerol, pentanols, n-octane, 2,2,4-trimethyl-pentane, lauryl alcohol 1-decanol, and silanols can be investigated. Each of these solvents is simulatable without expansion of the published OPLS-AA force field. Furthermore studies of mixed solvents in with polymer that include "dirt" materials such as butyric acid, propanoic acid, and 3-methyl butanoic acid would allow for the determination of dirt encapsulation effectiveness.
There are a large number of commercially relevant molecules which have yet to be simulated. The methods available in the GROMACS community combined with bash Unix scripting offers rapid, readily scalable approaches to deriving a large variety of their properties. Such results can be made available to the academic community to offer new insights in the correctness of force field parameterizations and eventually molecular design. These include ternary mixtures, block copolymers, and high molecular weight polymers.
Using the pressure tensor information and a code to determine the phase states of individual molecules, it will become possible to determine classical chemical engineering thermodynamic laws from computational molecular dynamics. These include Clausius-Clayperon, Antoine and Raoult's laws. After these are derivable via simulation, the psychometric chart could be reproduced.

Understanding reactivity from molecular dynamics: An insufficiently solved problem
One type of computational molecular dynamics that has not been explored in as much depth as could be expected. That is the idea of incorporating bond breaking and forming events into the classical ensemble. This could be accomplished by having a section of the force field control that checks if some criteria is met for an atom to transfer between molecules. Potential criteria include: (i) distances between atoms (ii) relative forces between atoms and (iii) bond distortion ratios relative to both the present and plausible future structure. The computational requirement for such a simulation would require that the cut-off distance for such checks be very short: in the 50-70 pm regimes. Although most bonds are approximately 15 pm, some are as large as 40 pm. In order for reactions to occur it is critical that all possible permutations of mixtures of topologies of the input molecules to be parameterized, otherwise the simulation would be highly unrealistic. Those molecules without parameterized structures would not have associated possible products from the liquid reaction. Often the dihedral interactions are the interaction which is not parameterized for a given bonded atom pairing combination within the OPLS-AA standard set. This boils down to a quantum chemistry problem.
The other consideration is what types of possible properties could be derived from such simulations. These would include a better understanding of the molecularity and concentration dependence of reaction rates from first principles, as well as the potential to determine the fundamental rate constant from a classical simulation. Thinking more broadly, another possible application is the determination of product distributions for various chemical reactions. One conceptual challenge that the author does not know how to approach is the determination of whether a single bond transitions to a bond break or a new intramolecular bond formed, or whether a new single bond, double or a triple bond is formed within a molecule or an adjacent one.

Conclusions
Firstly, some recent research was summarized and some modifications to it were suggested. Some of these expansions of work will allow new computational methods to become predictive of classical thermodynamic and kinetic quantities observed in experiments. Others will allow for chemical reactions in the liquid phase to take place. Although it is not clear of the OPLS-AA force field will offer accurate predictive capabilities for reaction product distributions, the key programmable logic checks for such software have been listed. While a few suggested problems flow naturally from recently published work, some are more difficult and overarching. These latter types of problems require researchers must combine their expertise to solve these newly suggested and difficult problems.