Computer simulation of disordered phases and interfaces
Projects:
1. Investigation of the intrinsic liquid-liquid and liquid-vapor interfaces
Intensive investigation of processes occurring at fluid interfaces, aiming at a molecular level understanding of the underlying phenomena came to the focus of scientific research in the past two decades, due to the development of experimental methods (e.g., nonlinear spectroscopy or reflectometry) that can distinguish between interfacial and bulk phase molecules, and to the rapid and continuous increase of the available computing capacities, which now enables computer simulation to meaningfully complement experimental investigations of fluid interfaces. In such simulations, however, one has to face the problem that, on the atomistic length scale, the liquid surface is corrugated by capillary waves, and hence finding the exact location of the interface is not a trivial task. We developed a method for identifying the truly interfacial molecules (ITIM), and hence locating the real, intrinsic surface of the liquid phase [1]. The method was later generalized also to non-planar surfaces [2]. The ITIM method was shown to be a nearly ideal compromise between computational cost and accuracy [3]. Having the intrinsic interface determined, profiles of various thermodynamic quantities (e.g., density, energy, free energy, pressure, surface tension, etc.) can be calculated relative to this intrinsic surface, and the thermodynamics of the surface can be characterized by these intrinsic profiles [4]. One of the important features of the ITIM method is that once the surface layer is fully determined, repeating the algorithm by neglecting the surface molecules will readily provide the list of molecules that form the second subsurface layer. This way, molecules constituting the subsequent subsurface layers can be detected, which allows an alternative way of describing surface thermodynamics through calculating the contribution of the subsequent layers to the thermodynamic quantity in question [4].
Investigation of near-critical fluid interfaces and fluid phases
Upon approaching the critical point (or, in the case of liquid-liquid interfaces, the critical mixing temperature), all properties of the two coexisting phases become increasingly similar to each other, while the surface tension approaches zero and the interfacial area increases. In the case of liquid-vapor interfaces, increasing amount of bubbles and liquid droplets occur in the liquid and vapor phase, respectively, and the densities of the two phases get closer to each other, while in the case of liquid-liquid interfaces the mutual solubilities of the two components increase, and hence the compositions of the two phases get more and more similar to each other. At the close vicinity of the critical point, these bubbles and liquid droplets become progressively less distinguishable from the increasing density fluctuations inside the bulk phases, and hence the detection of the interface as a molecularly sharp layer becomes more and more difficult as the concept of the interface itself gradually loses its physical meaning. In any case, the determination of the intrinsic surface in such systems requires a pre-sorting of the particles according to the phase they belong to. For this purpose, we proposed an algorithm based on finding the largest cluster of the like molecules [5], which was later further elaborated by my collaborators [6]. We try to push the limit of the detection of the interface as close to the critical point as possible, and develop a method that relates the determination of the intrinsic interface to the analysis of the bulk phase fluctuations in near-critical systems.
Relation between the intrinsic surface and the Gibbs dividing surface.
The use of the equimolar or Gibbs dividing surface is a well established concept in surface thermodynamics. Thus, in calculating the surface excess of any extensive quantity one has to define first a two-dimensional dividing surface between the two phases, and the resulting surface excess values depend on this choice. The surface excess values can be related to meaningful thermodynamic quantities through the particular choice of this surface as the Gibbs dividing surface, in which case, e.g., the surface excess Helmholtz free energy becomes identical with the surface tension, and the surface excess entropy with its temperature derivative. According to the general notion, only one particular choice of the dividing surface bears such property. We work on demonstrating that (i) there are infinite possible Gibbs dividing surfaces between two fluid phases, although all of them are of different shape, and hence there is only one such surface that is planar, a condition that has implicitly always been assumed but never explicitly stated, and (ii) demonstrate that the intrinsic covering surface of a fluid condensed phase, accessible by the ITIM method or its generalized version GITIM, is always a Gibbs dividing surface.
Single particle dynamics at near-critical interfaces
The identification of the full set of the molecules that constitute the liquid surface at every instant opens also the possibility of studying the dynamical properties of the interfacial molecules. We studied single particle dynamics both at liquid-vapor [7,8] and liquid-liquid [9] interfaces in systems where practically no particle penetrated to the opposite phase in the simulations. Now we would like to extend this study to near-critical interfaces, and address the conceptual questions arising from the increasing similarity of the composition of the two phases.
Collective dynamics at the free liquid surface.
We try to analyze the collective dynamics of the molecules at neat liquid surfaces. Unlike in the case of single particle dynamics, here one does not have to face to the problem of time scales, since accessing the propagation of different waves requires the knowledge of the identity, position and speed of the full set of particles at every instant, but this set does not have to be the same from instant to instant. We are interested in how the longitudinal and transverse currents and their spectra differ at the liquid surface from those in the bulk liquid phase, in the conceptual issues stemming from the fact that currents along the macroscopic surface plane and its normal axis have to be distinguished, and in the question whether there is any meaningful, detectable surface current along the surface normal axis.
Determination of the surface tension contribution of various moieties at fluid interfaces
At the molecular level, surface tension originates from the fact that individual particles have to pay energetic and entropic penalties for being at the surface of their phase. This fact naturally gives rise to the question of how the different molecules contribute to the interfacial tension. Since the surface tension itself is related to the imbalance between the lateral and normal components of the pressure along the interface normal, and the latter component has to be, due to the requirement of the mechanical stability, constant, the surface tension of the system can be accessed through the lateral pressure profile. We proposed a computationally very efficient and accurate way of calculating the lateral pressure profile by distributing the lateral pressure contributions among the individual atoms [10]. This treatment opened the possibility of calculating the surface tension contribution of various molecules and moieties [11]. With the aid of this method, we try to address the following issues.
(1) Applying this method for the free surface of aqueous solutions of surfactants [12] we got the unexpected result that surfactant counterions have a huge surface tension contribution, this contribution is even of different sign for Na+ (next to the dodecylsulfate (DS-) anion) and for Cl- (next to the dodecyl-trimethylammonium (DTA+) cation), and the contributions of the surfactant heads and counterions are of opposite sign. This finding leaves the floor open to several explanations. Thus, the surface tension contribution of the ions might be related to their softness/hardness (according to the Hofmeister series), but also to their net charge (through the geometric asymmetry of the charge distribution of the water molecules), an effect that was shown to have important consequences even in the surface affinity of alkali and halide ions [13].
(2) Although interfacial tension arises from the fact that at least one of the two coexisting phases experiences, with respect to the bulk liquid phase, a less favourable environment near the interface, it might well be the case that the other component is in a more favourable environment at the vicinity of the opposite phase than in its bulk phase. Interfaces between water and various apolar liquids might be possible examples for such systems, because the apolar molecules are not expected to experience an energetically much different environment at the aqueous interface than in their own bulk liquid phase. We investigate how the molecules forming the two opposite phases contribute to the interfacial tension at liquid-liquid interfaces between water and various apolar organic liquids.
2. Investigation of the molecular mechanism of anesthesia
Although the effect of anesthesia is used in surgery for more than a century, and it is also long known that anesthetics are assembled in the membrane of the cells, its molecular mechanism is still largely unknown. It should also be noted that lipid membranes have two biologically relevant phases: the liquid crystalline () and gel () phases, the former existing at lower pressures and higher temperatures. According to long standing experimental observations, the effect of anesthesia is reverted at high (100-1000 bar) pressures, therefore, any explanation of its molecular mechanism should also account for the pressure reversal. It is also worth noting that general anesthetics show a great chemical variety, and hence no specific interaction can account for anesthetic effects. Existing conjectures concerning the molecular level explanation of the phenomenon of anesthesia can be divided to two major groups. Lipid theories assume that anesthesia is caused by the fact that anesthetic molecules perturb certain properties of the cell membrane, whereas protein theories assume that the changes in the membrane structure, caused by the anesthetics, alter the conformation of certain transmembrane proteins. However, the large chemical variety of general anesthetics practically excludes any general mechanism, valid for all anesthetics, in the frame of the protein theories. Further, it is very difficult to imagine any evolutionary mechanism that led to the development of any specific interaction of proteins with a chemically totally inactive species, such as xenon, which is not even present in a considerable concentration on Earth.
In the frame of the lipid theories, several membrane properties (e.g., lateral density, membrane fluidity) have been suggested to lay behind the mechanism of general anesthesia, however, these properties have not been unambiguously identified yet. In 1997, Cantor presented thermodynamic arguments demonstrating that if anesthesia is caused by the switch between two conformations of a membrane-bound protein, and anesthetics alter the lateral pressure profile across the membrane in such a way that it changes the cross-section area profile of both conformers of this protein non-uniformly, the ratio of the active and passive conformers depends exponentially on the anesthetic concentration [14]. This finding made the lateral pressure profile to be a prominent candidate of being such a property. As an alternative explanation, Heimburg presented a consistent thermodynamic formalism [15], and assuming that anesthetics dissolve ideally in the , but do not dissolve at all in the phase, they showed that the observed temperature shift of the gel-liquid crystalline phase transition is analogous with the well known freezing point depression of simple solutions. Since freezing point depression is a colligative property, its magnitude solely depends on the number of dissolved particles, being independent from their chemical nature. Thus, given that equal amount of different anesthetics are used, the shift of the - phase transition temperature they induce depends only on their partition between the aqueous phase and the membrane, in accordance with well known experimental facts. This way, pressure reversal can easily be explained, as pressure shifts the gel-liquid crystalline phase transition to higher temperatures. This theory is also able of explaining several subtle points of anesthesia, such as the cut-off effect of n-alkanols.
We try to identify such properties using the following criteria: (i) all the chemically different anesthetic molecules considered should alter the property in question in the same way, and (ii) the increase of the pressure should alter the property in question in the opposite way. We found [16,17] that anesthetics prefer, among others, to stay in the outer part of the apolar hydrocarbon region, close to the dense region of the polar headgroups. The anesthetic molecules located in this region induce a lateral swelling of the membrane by pushing the lipid chains farther away from each other, and thus creating some additional empty space here. This additional space causes (i) a natural increase of the free volume fraction in this region, (ii) increases the lateral mobility of the lipid molecules, and (iii) decreases the lateral pressure in the nearby region of the amide and ester groups, i.e., at the boundary of the apolar and polar parts of the membrane, where the anesthetic molecules already cannot penetrate, but their vicinity still pushes, on average, farther away the lipid molecules from each other in this crowded region. Further, since the lateral swelling of the membrane is clearly reverted by the increase of the pressure, all these changes are also pressure reversible. However, although these results make the lateral density a strong candidate for being the relevant membrane property (and being the cause of the changes in the lateral pressure profile), this hypothesis (as any other one) needs to be tested also by simulating molecules that are chemically similar to known anesthetics but do not bear anesthetic properties, and demonstrating that, unlike anesthetics, they do not induce such changes in the membrane structure. In this field, we try to address the following issues.
(1) Although we have checked in our preliminary study several, chemically rather different anesthetics (e.g., halothane, enflurane, chloroform, diethyl ether, sevoflurane), the effect of probably the most mysterious anesthetic, namely xenon, on the membrane properties also needs to be tested. Xe, being a noble gas, is not supposed to have much interaction beyond the hard core one with the molecules of the membrane. The behavior of Xe in lipid membranes should thus be a stringent test of our above hypothesis as well as of those of Cantor and Heimburg.
(2) We test our above hypothesis by replacing anesthetics with chemically similar non-anesthetics (e.g., diethyl ether with pentane, chloroform with CCl4, Xe with Ne, etc.) to see whether such non-anesthetic molecules also induce lateral swelling of the membrane, leading to all other related changes, or not.
(3) We test the hypothesis of Heimburg by calculating the lateral density in membrane simulations without embedded molecules, with anesthetics, and also with non-anesthetics at several temperatures, both below and above the -phase transition, and estimate the temperature of the phase transition, in order to see whether anesthetics indeed shift the - phase transition temperature to lower values, and also check whether no such shift is seen in the case of non-anesthetics.
3. Investigation of problems related to solid interfaces
Coating of magnetite nanoparticles in aqueous environment
Dispersing magnetite nanoparticles in a fluid phase is a conventional way to prepare magnetic liquids. Such magnetic liquids are not only of fundamental physico-chemical interest, but also of great practical importance. Besides several important technological applications, the potential use of magnetic liquids in medical applications is also the subject of intensive scientific research. Thus, by applying external magnetic field, drug molecules that are attached to the surface of the magnetite nanoparticles can be selectively brought to the targeted area of the body. However, the dispersion of magnetite has to be kinetically stabilized, and the surface of the nanoparticles has first to be activated by appropriate non-toxic, physiologically neutral surfactants for the attachment of the selected drug molecules. For designing such coatings, the role of the various intermolecular interactions, acting at the surface of the magnetite nanoparticle, has to be deeply understood first. For this purpose, computer simulation seems to be a particularly suitable tool. As a first step in this direction, we have already calculated the adsorption isotherm of water on magnetite, which proved to be in an excellent agreement with experimental data [18]. Although this result gives us, in general, confidence in the molecular model used, this model does not take into account chemisorbed water molecules at the magnetite surface, a key factor when studying surface-adsorbate interactions. For this purpose, we work on determining the number and position of these chemisorbed water molecules at the surface of the unit cell of the magnetite crystal by extensive ab initio calculations. Having the chemisorbed water molecules accurately located, fractional charges can be assigned to them, and the previously successfully used magnetite model can be be completed by these chemisorbed surface waters. The new magnetite surface created this way needs to be be tested by simulating its aqueous surface, calculate various properties of the surface water molecules, and compare the results with existing experimental data. After verifying the model, we are going to investigate the adsorption of several possible coating agents of various chemical nature at the surface of magnetite from the aqueous phase. Among these coating agents, malonic acid and citric acid are small molecules with several functional groups, which are expected to bind to the magnetit surface by forming several H-bonds, oleic acid is an amphiphilic molecule forming a bilayer at the magnetit surface, whereas polyethylene glycol is a conventionally used amphiphilic polymer.
Adsorption of various small molecules on amorphous ice under interstellar conditions
The interaction of small molecules in the interstellar medium (ISM) with the surface of amorphous (LDA) ice is an important issue in astrochemistry, because the surface catalyzed reactivity of these molecules play a non-negligible role in the chemistry of the ISM. Further, several of these reactions might lead to the formation of precursors of large biomolecules, and hence could have possibly played a role in the prebiotic evolution. Experimental identification of small molecules in the ISM is a rapidly developing area, making computational modeling of the behavior of these compounds under astrophysical conditions also an increasingly important field of chemical research. While the reactivity of these molecules under such conditions is the subject of intensive investigations by quantum chemical methods, the basic question whether these molecules do adsorb on LDA ice under interstellar conditions, and hence such quantum calculations correspond to relevant arrangements of the molecules is much less studied. We study the adsorption of a number of small molecules, existing in the ISM, at the surface of LDA ice. These molecules include benzonitrile, a molecule bearing a large (i.e., about 4D) dipole moment, acetamide, a H-bond donor molecule, cyanamide, which can act both as a H-donor and a H-acceptor, and propylene oxide, the only chiral molecule that has been identified in the ISM so far.
4. Investigation of problems related to bulk fluid phases
Investigation of the thermodynamics of mixing of various compounds
We developed a method to calculate the Helmholtz free energy, energy, and entropy of mixing of two compounds in computer simulation [19] using the method of thermodynamic integration [20]. This method has mostly been used so far to test the performance, including the miscibility, of possible model combinations against existing experimental data. However, having, besides the free energy, also the energy and entropy of mixing determined in the entire composition range this way, the analysis of the energetic and entropic contributions can provide a detailed insight into the thermodynamic background of the miscibility (or non-miscibility) of the two compounds, which is clearly related to the local structure of the constituting molecules [21]. We are applying this method to study the mechanism of miscibility of various H-bonding, aprotic polar and apolar compounds (e.g., water with methanol and with DMF or acetonitrile, CCl4 with methanol, acetonitrile and hexane, hexane with hexanol) to get a deeper understanding of the thermodynamic background of the mixing properties of compounds of different chemical characters.
Investigation of potential green solvents
In the present days, there is a considerable effort of substituting conventional, toxic solvents by environmentally more friendly alternatives in industrial processes. Such solvents include room temperature ionic liquids (RTILs), mixtures of supercritical CO2 (scCO2) and polar co-solvents (e.g., acetone, methanol), and deep eutectic solvents. The properties of these solvents can be fine-tuned due to the delicate interplay of the polar or charged and apolar groups of the constituting molecules, which results in the formation of various microscopic polar and apolar aggregates of the constituting particles or moieties. We are studying the local structure in ionic liquids consisting of an alkyl-3-methylimidazolium cation and a quasi-spherical anion, and in their mixtures with polar co-solvents, in mixtures of scCO2 with methanol and acetone in the entire composition range, and in the deep eutectic mixtures of choline chloride and urea, as a function of the composition and, in the case of RTILs, also of the alkyl chain length by computer simulation. Self-association and local structure are characterized by determining cluster size distributions and Voronoi analysis, and the results are related to the known solvation properties of these systems.
References
[1] Pártay, L. B.; Hantal, Gy.; Jedlovszky, P.; Vincze, Á.; Horvai, G. J. Comp. Chem. 2008, 29, 945.
[2] Sega, M.; Kantorovich, S.; Jedlovszky, P; Jorge, M. J. Chem. Phys. 2013, 138, 044110.
[3] Jorge, M.; Jedlovszky, P.; Cordeiro, M. N. D. S. J. Phys. Chem. C. 2010, 114, 11169.
[4] Sega, M., Fábián, B., Jedlovszky, P. J. Chem. Phys. 2015, 143, 114709.
[5] Pártay, L. B.; Horvai, G.; Jedlovszky, P., J. Phys. Chem. C. 2010, 114, 21681.
[6] Sega, M.; Hantal, Gy. Phys. Chem. Chem. Phys. 2017, 19, 18968.
[7] Fábián, B.; Senćanski, M. V.; Cvijetić, I. N.; Jedlovszky, P.; Horvai, G. J. Phys. Chem. C. 2016, 120, 8578.
[8] Fábián, B.; Horvai, G.; Sega, M.; Jedlovszky, P. J. Phys. Chem. B. 2017, 121, 5582.
[9] Fábián, B.; Sega, M.; Horvai, G.; Jedlovszky, P. J. Phys. Chem. C. 2020, 124, 2039.
[10] Sega, M. Fábián, B.; Jedlovszky, P. J. Chem. Theory Comput. 2016, 12, 4509.
[11] Sega, M. Fábián, B.; Horvai, G.; Jedlovszky, P. J. Phys. Chem. C. 2016, 120, 27468.
[12] Hantal, Gy.; Sega, M.; Horvai, G.; Jedlovszky, P. J. Phys. Chem. C. 2019, 123, 16660.
[13] Hantal, Gy.; Horváth, R. A.; Kolafa, J.; Sega, M.; Jedlovszky, P. J. Phys. Chem. B. 2020, 124, 9884.
[14] Cantor R. S. Biochemistry 1997, 36, 2339.
[15] Heimburg, T.; Jackson, A. D. Biophys. J. 2007, 92, 3159.
[16] Fábián, B.; Sega, M.; Voloshin, V. P.; Medvedev, N. N.; Jedlovszky, P. J. Phys. Chem. B 2017, 121, 2814.
[17] Hantal, Gy.; Fábián, B.; Sega, M.; Jójárt, B.; Jedlovszky, P. Biochim. Biophys. Acta 2019, 1861, 594.
[18] Tombácz, E.; Hajdú, A.; Illés, E.; László, K.; Garberoglio, G.; Jedlovszky, P. Langmuir 2009, 25, 13007.
[19] Darvas, M.; Jedlovszky, P.; Jancsó, G. J. Phys. Chem. B 2009, 113, 7615.
[20] Mezei, M.; Beveridge, D. L. Ann. Acad. Sci. N.Y. 1986, 482, 1.
[21] Idrissi, A.; Polok, K.; Barj, M.; Marekha, B.; Kiselev, M.; Jedlovszky, P. J. Phys. Chem. B 2013, 117, 16157.