Estrategias de muestreo en simulación molecular con cuerpos rígidos : evaluación crítica y aplicación al cálculo de energía libre de solvatación

Abstract

La presente investigación se centra en dos ejes principales. El primero está relacionado con metodologías de muestreo para la simulación de moléculas modeladas como cuerpos rígidos. Utilizando dichas metodologías fue posible desarrollar un segundo eje de investigación, el cual tiene como objeto de estudio a los cálculos de energía libre de solvatación y de propiedades termodinámicas relacionadas. El Capítulo 3 tiene como objeto de estudio al método de Monte Carlo Híbrido (MCH) aplicado a la simulación de cuerpos rígidos. Se realizó un estudio sistemático de la eficiencia de MCH, y se desarrolló un procedimiento de optimización en línea que permite obtener parámetros óptimos o cuasi-óptimos de MCH para los sistemas en fase líquida investigados: propano, benceno y agua. La implementación de MCH, realizada en LAMMPS, fue motivada por aspectos objetables del método del ensamble expandido: usualmente se combinan pasos de DM y MC violando la condición del balance detallado de probabilidad y, además, no se realiza una equilibración luego de las transiciones entre los macroestados. En el Capítulo 4 se realiza una validación empírica del método del ensamble expandido, para lo cual se determinó un perfil de energía libre “correcto”, que se obtuvo a partir de muestras independientes de cada macroestado empleando MCH. Asimismo, se evaluaron diferentes frecuencias de transición entre los macroestados, y en todos los casos el perfil de energía libre coincide con el resultado de referencia. De este análisis se desprende que los aspectos objetables no introducen errores sistemáticos importantes, y se concluye que el muestreo secuencial del ensamble expandido con DM y MC arroja resultados confiables. La validación del método del ensamble expandido permitió abordar el desafío de predecir la solvatación en líquidos iónicos y calcular propiedades macroscópicas requeridas en el diseño de procesos químicos, tal como se detalla en el Capítulo 5. El CO2 proveniente de procesos de combustión constituye una de las fuentes principales de gases de efecto invernadero de generación antropogénica y el potencial de [emim][B(CN)4] para su captura ha sido recientemente reportado. Con fines comparativos, se analizaron los solventes [emim][B(CN)4] y [emim][Tf2N], tanto en la condición de dilución infinita, como en altas concentraciones. Las simulaciones confirmaron lo observado experimentalmente, en tanto se requiere de una presión considerablemente menor en [emim][B(CN)4] para obtener una misma molalidad de CO2. Por otra parte, se predijeron los coeficientes de actividad a dilución infinita de los solutos hexano, benceno, ciclohexano, etanol y agua en [emim][B(CN)4], y se compararon con valores experimentales reportados recientemente en la literatura. Finalmente, en el Capítulo 6 se exponen las conclusiones y consideraciones finales de la presente investigación, como también algunas sugerencias para trabajos futuros. En el Capítulo 1 se introducen conceptos de Mecánica Estadística y Simulación Molecular (SM) que han sido fundamentales para el desarrollo de esta tesis. A esto antecede una revisión bibliográfica en la cual se destacan factores determinantes para que la Simulación Molecular constituya actualmente una herramienta de investigación en la Ingeniería Química. En el Capítulo 2 se analizan y desarrollan métodos de Dinámica Molecular (DM) de cuerpos rígidos en los ensambles NVE y NVT, motivado principalmente por resultados insatisfactorios obtenidos con un integrador NVT implementado en softwares muy utilizados por la comunidad de SM, y que acopla termostatos independientes a los grados de libertad traslacionales y rotacionales. En particular, para el sistema de agua líquida, dicha estrategia no logra reproducir la temperatura especificada, inclusive para un paso de tiempo pequeño de 1 fs. Ante esto, es natural preguntarse si el problema podría estar vinculado con el integrador NVE, que es la base del método NVT. Una anomalía detectada fue la violación del teorema de equipartición de la energía cinética, lo cual tiene origen en la discretización de la integración numérica. Esto ya ha sido reportado en la literatura, pero sólo se mostró el efecto de esos errores sobre diversas propiedades termodinámicas y, hasta donde se sabe, no existen trabajos que hayan revisado el método NVT mencionado. Esto motivó la derivación de un integrador NVE con una formulación matemática simplificada de las ecuaciones de movimiento. En la presente investigación se demuestra que el origen de las anomalías está relacionado con la disparidad de los momentos de inercia. Además, se presenta un procedimiento simple para disminuir el efecto de los errores de discretización, que viene dado por la redistribución de las masas atómicas.

The present investigation focuses on two main topics. The first is related to sampling strategies for the simulation of molecules modeled as rigid bodies. This trigger a second line of research, which aims to calculate solvation free energies and related thermodynamic properties. Chapter 1 introduces concepts of Statistical Mechanics and Molecular Simulation (MS) needed for the development of this thesis. This is preceded by a literature review that highlights the major factors that drove Molecular Simulation to become a research tool in Chemical Engineering. Chapter 2 is devoted to the analysis and development of integrators for Molecular Dynamics (MD) with rigid bodies in the NVE and NVT ensembles. This was mainly motivated by unsatisfactory results obtained with an NVT integrator that couples independent thermostats to the translational and rotational degrees of freedom, which is implemented in software packages widely used by the MS community. In the case of liquid water, such strategy fails to reproduce the specified temperature, even for a time step as small as 1 fs. Thus, it is natural to ask whether the problem could be linked to the NVE integrator, which is the basis of the NVT method. An anomalous violation of the kinetic-energy equipartition theorem was detected, originated by the time discretization employed in the numerical integration. This has already been reported in the literature, but the discussion was limited to the effect of these errors on various thermodynamic properties and, as far as is known, no one has reviewed the NVT method mentioned above. This motivated the derivation of an NVE integrator with a simplified mathematical formulation of the equations of motion. In the present investigation it is shown that the origin of the anomalies is related to a large difference in the moments of inertia of a rigid molecule. In addition, a simple procedure for reducing the effect of discretization errors is presented, which consists in redistributing the atomic masses. In Chapter 3, the Hybrid Monte Carlo (HMC) method is applied to the simulation of rigid bodies. A systematic study on the efficiency of HMC was carried out. Then, an on-the-fly optimization procedure was developed, which allows obtaining optimal or quasi-optimal parameters for HMC. The liquid-phase systems assessed were propane, benzene, and water. The implementation of HMC, carried out in LAMMPS, was motivated by questionable aspects of the expanded ensemble method: MD and MC steps are usually combined in a way that violates the condition of detailed balance of probabilities and, in addition, equilibration is not performed after a transition between macrostates takes place. In Chapter 4, an empirical validation of the expanded ensemble method is shown, for which a “correct” free energy profile was determined from independent samples drawn at each macrostate using HMC. Likewise, different transition frequencies between macrostates were evaluated, and in all cases the free energy profile agrees with the reference result. From this analysis, it seems that the questionable aspects do not introduce important systematic errors and the sequential sampling of an expanded ensemble alternating MD and MC yields reliable results. The validation of the expanded ensemble method allowed to tackle, in Chapter 5, the challenge of predicting solvation in ionic liquids and macroscopic properties required in the design of chemical processes. Carbon dioxide produced in combustion processes is one of the main sources of anthropogenic generation of greenhouse gases. On the other hand, the potential of [emim][B(CN)4] for its capture has recently been reported. For comparison, the solvents [emim][B(CN)4] and [emim][Tf2N] were analyzed, both in the condition of infinite dilution and at high concentrations. The simulations confirmed what was observed experimentally, that is, a considerably lower pressure is required in [emim][B(CN)4] to obtain the same molality of CO2. The infinite-dilution activity coefficients of solutes hexane, benzene, cyclohexane, ethanol, and water in [emim][B(CN)4] were also calculated, and a comparison was made with experimental values recently reported in the literature. Finally, Chapter 6 presents the conclusions and final considerations of the present investigation, as well as some suggestions for future work.