STATISTICAL M E C H A N I C S OF SOLVATION OF M A C R O PARTICLES By Reka Z. Vasarhelyi B . Sc. (Chemistry) University of British C o l u m b i a , 1990 A THESIS THE S U B M I T T E D IN P A R T I A L F U L F I L L M E N T REQUIREMENTS FOR T H E DEGREE MASTER OF OF SCIENCE in THE FACULTY O F G R A D U A T E STUDIES CHEMISTRY We accept this thesis as conforming to the required standard THE U N I V E R S I T Y O F BRITISH COLUMBIA September 1992 © Reka Z. Vasarhelyi, 1992 OF In presenting this degree at the thesis in University of partial of this department publication or and his or her permission. Chemistry The University of British Columbia Vancouver, Canada Date DE-6 (2/88) purposes may representatives. of this thesis for financial gain Department of Oct. 13 1992 the requirements study. I further agree that thesis for scholarly by of British Columbia, I agree that the freely available for reference copying fulfilment be It shall not be advanced Library shall make it permission for extensive granted is for an by understood allowed the that head of my copying or without my written Abstract Integral equation theories are employed to study the solvation of large particles. The solutions studied consist of spherically symmetric solvent particles which interact through hard sphere and Lennard-Jones potentials, and hard sphere solutes. In particular, the Ornstein-Zernike (O-Z) equation is solved with the hypernetted-chain ( H N C ) closure, to obtain the pair correlation functions of mixtures at infinite dilution. T h e pair correlation functions in the O-Z equation and the H N C closure are expressed as power series in solute density to yield a pair of coupled equations which determine the derivatives with respect to solute density of the solvent-solvent pair correlation functions. T h e latter describe the perturbation of the solvent upon addition of a single solute particle. T h e derivatives are analysed to yield components that scale as the volume and surface area of the macroparticle, and which are then identified as changes i n solvent structure due to the presence of a finite size particle and a flat surface respectively. F r o m the pair correlation functions and their derivatives the excess internal energy, Helmholtz free energy, and entropy of solvation are calculated. T h e excess quantities are also separated into contributions from finite size and surface effects. B o t h components of the excess internal energy are negative at low densities, and become less negative for high density liquids. T h e magnitude and sign of the two contributions to the energy depend on physical conditions such as temperature and pressure. T h e excess entropy of solvation is negative for all systems studied, indicating that introduction of a macroparticle or a flat surface increases the spatial ordering of a bulk liquid. n Table of Contents Abstract ii List of Tables v List of Figures vi Acknowledgement viii 1 Introduction 1 2 Theory 8 2.1 The R a d i a l Distribution Function 8 2.2 Choice of M o d e l 9 2.3 Thermodynamics 11 2.4 Integral equations 15 2.4.1 The Ornstein-Zernike equation 15 2.4.2 Choice of Closure 18 2.4.3 Numerical Details 22 3 4 Results and Discussion 23 3.1 Input Parameters 23 3.2 Radial Distribution Functions 24 3.3 Thermodynamics 42 Conclusions 53 iii Appendices 56 A 56 H e l m h o l t z Free E n e r g y for a M i x t u r e Bibliography ^QQ iv List of Tables 3.1 Input Parameters; H S / H S refers to a hard sphere solvent with a hard sphere solute; I I S / L J refers to a Lennard-.Jones solvent w i t h a hard sphere solute 3.2 23 Excess Energy of Solvation from extrapolation of < U ex > /N^d 2 vs. d plots 47 3.3 Excess Energy of Solvation from integration of <*)</!(r)and 8g2(r) 47 3.4 Excess Energy of Solvation at constant total density from extrapolation of < U ex 3.5 > /N^P vs. d plots 49 Excess Energy of Solvation at constant total density, from integration of (fyl(r)and 8g2{r) 49 3.6 Excess Free Energy of Solvation, from E q . 2.26 50 3.7 Excess Free Energy of Solvation, from E q . 2.61 50 3.8 Excess Entropy of Solvation 52 List of Figures 3.1 Solute-solvent radial distribution functions for H S / H S solutions; p* = 0.8. T h e solid, dotted and dashed curves correspond to solvent-solute diameter ratios of 1:1, 1:10, and 1:100, respectively 3.2 25 Solute-solvent radial distribution functions for H S / L J solutions; p* = 0.6, T* = 1.35. T h e solid, dotted, and dashed curves correspond to solventsolute diameter ratios of 1:1, 1:10, and 1:100, respectively 3.3 Solvent-solvent radial distribution function for H S / L J solution; p* = 0.6, T* = 1.35 3.4 28 Derivatives of solvent-solvent radial distribution functions for H S / H S solutions; p* = 0.8. T h e solid, dotted, and dashed curves correspond to solvent-solute diameter ratios of 1:1, 1:10, and 1:100, respectively 3.5 26 29 Derivatives of solvent-solvent radial distribution functions for H S / L J solutions; p* = 0.6, T* = 1.35. The solid, dotted, and dashed curves correspond to solvent-solute diameter ratios of 1:1, 1:10, and 1:100, respectively. 30 3.6 8g (r) vv at r = Id for the H S / L J solution; p* = 0.6, T* = 1.35 31 3.7 Sg^(r) at r = 5rf for the H S / L J solution; p" = 0.6, T* = 1.35 32 3.8 Sg„ (r) at r = lOd for the H S / L J solution; p* = 0.6, T* = 1.35 33 3.9 6g {r) at r = 30rf for the H S / L J solution; p" = 0.6, T* = 1.35 34 v vv 3.10 Sgl(r) for H S / L J solution; p* = 0.6, T* = 1.35. T h e solid, dotted, and dashed lines correspond to solutions of Eqs. 3.61 using different solute diameter pairs 36 vi 3.11 8g2(r) for H S / L J solution; p' = 0.6, T" = 1.35. T h e solid, dotted, and dashed lines correspond to solutions of Eqs. 3.61 using different solute diameter pairs 3.12 6gl(r) 37 for H S / L J solution on an expanded scale; p* = 0.6, T* = 1.35. . . 3.13 8g2(r) for H S / L J solution on an expanded scale; p* = 0.6, T* = 1.35 . . . 3.14 8g (r) vv 38 39 at constant total density for H S / L J solution; p = 0.6, T* = 1.35. T h e solid, dashed, and dotted curves correspond to solvent-solute diameter ratios of 1:1, 1:10, and 1:100, respectively 40 3.15 Excess Energy of Solvation as a Function of Solute Diameter 43 3.16 Excess Free Energy of Solvation as a Function of Solute Diameter 44 3.17 Excess Entropy of Solvation as a Function of Solute Diameter 45 VI1 Acknowledgement I would like to acknowledge the help and support I have received from several people during the course of this work. M y supervisor D r . G . N . Patey suggested an interesting topic, and provided expert guidance throughout. My fellow workers contributed to an excellent work environment under the palm tree. Discussions with D r . P h i l i p A t t a r d and Daniel Berard were especially helpful and enjoyable. Finally, I would like to thank my husband P a u l , who was, as always, understanding and supportive of my efforts. v 111 Chapter 1 Introduction A great deal of recent work in liquid state theory has focussed on inhomogeneous liquids. Such systems include dense fluids in contact with a wall or macroparticle, and liquids i n slits. T h e elucidation of the microscopic structure of a liquid in contact w i t h a macroscopic, body is important for the understanding of electrode processes, the phenomenon of wetting, and interactions between large particles. T h e general purpose of any statistical mechanical theory is the development of microscopic models which reproduce macroscopic properties such as thermodynamic quantities. A theoretical treatment of the liquid state is unique in the sense that a simple reference system, such as the harmonic solid and ideal gets, for solids and gases respectively, does not exist. Unlike a solid or gas, a real liquid cannot be thought of as a departure from an ideal state. In this thesis, only classical fluids are considered. T h i s means that the variables of kinetic energy, velocity, and momentum, are decoupled from those of potential energy, position, and orientation. Since the kinetic contribution to thermodynamic quantities is known, only the contribution of the potential energy will be studied here. Also, we restrict ourselves to simple liquids. Specilically, we mean dense fluids of monatomic, spherical, argon-like particles. Central to the microscopic description of liquids of simple spherical particles is the radial distribution function, g(r). T h i s function, which represents the probability of finding a second particle at a distance r from a particle at the origin, is particularly useful since it can be related to the structure factor. T h e latter is a physical observable 1 Chapter 1. 2 Introduction of X - r a y and neutron scattering experiments. Furthermore, if we make the simplifying assumption that the total potential energy of the system is a sum of pairwise additive potentials, then all thermodynamic quantities can be expressed as integrals over the radial distribution function. Examples of two such relations are the compressibility equation, (1.1) and the pressure equation, (1.2) which is also known as the virial equation of state [1]. If an exact g(r) were available, both equations would give the same value for the pressure. For approximate radial distribution functions, the discrepancy between the pressures calculated from the two equations is an indication of the thermodynamic consistency of the approximate theory [1]. T h e function g(r) can be calculated using two different, but complementary, techniques, computer simulations and analytical theories [2]. T w o types of computer experiments are generally performed for liquids. A Molecular Dynamics ( M D ) simulation involves the solution of the classical equations of motion for a set of particles, assigned initial velocities in accordance with the M a x w e l l - B o l t z m a n n distribution for a given temperature. T h i s technique is appropriate for the study of time-dependent processes. T h e Monte Carlo ( M C ) method takes a probabilistic, rather than a deterministic approach. F r o m an initial set of randomly assigned coordinates, further configurations are generated by random displacements. T h e configuration is accepted or rejected depending on the B o l t z m a n n factor, and therefore the probability of such an arrangement of the particles of a system at equilibrium. A M C simulation yields e q u i l i b r i u m ensemble averages. Since particles in a dense fluid are influenced by many other particles, the interaction between two molecules is a sum of pair and many-body forces. For argon-like liquids, potentials Chapter 1. Introduction 3 which accurately fit experimental data have been developed by Barker [3], using quantum mechanical perturbation methods. In most cases however, simple analytic expressions for the effective interaction between two particles are employed. A pair potential is assumed which retains the main features of the interparticle interaction. For simple liquids, the dominant features are a short-range strong repulsion and a long-range weak attraction. T h e potential energy of the model fluid is then taken to be a sum of pair interactions over a l l pairs of particles. In a computer simulation, the form of the pair-potential is part of the input, and is therefore defined exactly. That is, the results of a computer simulation are essentially exact for the model liquid. A n a l y t i c a l theories, on the other hand, involve further approximations. W h e n comparing the results of analytical calculations to experimental data, the source of deviations is ambiguous. In contrast, comparison with computer simulation tests the mathematical approximations of the theory. Current analytical theories can be placed in two broad categories, perturbation and integral equation theories. T h e former, due originally to the work of Zwanzig [4], treats a fluid whose molecules interact v i a a potential containing both a repulsion and an attraction as a perturbation of a fluid with a purely repulsive pair potential. W h i l e a reference system consisting of hard spheres is not an exact ideal state, it has proven to be a practical reference state. This is due to the fact that the structure of a fluid is determined almost entirely by the repulsive part of the potential. Also, hard sphere liquids have been completely characterised by computer simulations. Perturbation theories are easy to work out, but their success has been limited to simple systems [5]. Integral equations are somewhat more demanding, but have been more successful for complex systems. Integral equation theories are based on approximate equations for the pair distribution function. One of the first such relations was developed by K i r k w o o d in the 1930's [6]. T h e differentiation of the equilibrium n-particle density w i t h respect to a coupling parameter gives the K i r k w o o d hierarchical equations. If the differentiation is Chapter 1. Introduction 4 carried out w i t h respect to the coordinates of a given particle, another set of equations, called the B o r n - G r e e n - Y v o n hierarchy is obtained. As the names suggest, both theories consist of sets of coupled equations which give the n-particle distribution functions in terms of the (n+l)-particle distribution functions. Another set of integral equtions, which was derived in the 1950's, is based on a simple physical idea, proposed by Ornstein and Zernike several decades earlier [7]. T h e mathematical statement of their idea is, (1.3) This equation, and its relation to g(r), will be discussed in greater detail in the next section. For now, it is sufficient to say that the total influence of particles 1 and 2 on each other, given by the total correlation function / / ( r i , ! ^ ) , is a sum of a direct part described by the direct correlation function c ( r i , r 2 ) , and an indirect part described by the integral term in E q . 1.3. T h e quantity p represents the bulk density. For a homogeneous isotropic fluid, the interparticle forces are central; therefore, the correlation functions depend only on the magnitude of the separation between particles [7]. In such a case, the Ornstein-Zernike (O-Z) equation contains only two unknowns, namely h(r) and c(r). Hence, given another relation between the total and direct correlation functions, a complete solution can be obtained. A p p r o x i m a t e equations relating h(r) and c(r), such as the Perçus-Yevick ( P Y ) , mean spherical ( M S A ) , hypernetted- chain ( H N C ) , and reference hypernetted-chain ( R H N C ) approximations, combined w i t h the O - Z equations form several integral equation theories for g(r)T h e homogeneous integral equation theories have been extended to the case of a macroparticle in contact with a bulk fluid, by considering a l i m i t i n g case of a homogeneous two component mixture [8]. In the limit of large solute radius and infinite solute dilution, the mixture becomes a single component fluid in contact w i t h a flat surface [9]. Chapter 1. Introduction 5 Perram and W h i t e [8], have shown that it is physically reasonable to take this l i m i t by demonstrating that for a growing hard sphere solute, g(r) becomes v i r t u a l l y independent of the size of the macroparticle for solute sizes as small as fifteen times the solvent diameter. In a concurrent but independent study Henderson et cd. [9] evaluated analytically the infinite radius limit of the 0 - Z equation. They combined this result w i t h the P - Y approximation to obtain a one dimensional wall-solvent distribution function. T h e wall-solvent 0 - Z equation, in combination w i t h various closure relations, has been used to calculate the density profile of systems such as hard sphere and LennardJones fluids in contact with hard or attractive walls [10]. T h e technique has also been successfully applied to the problem of developing a molecular description of the electrical double layer [11, 12]. For liquids near a phase transition, density profiles calculated from the homogeneous 0 - Z equation do not agree well with those obtained from simulations [3], T h e inhomogeneous O-Z equation has been applied to such cases [13]. T h e inhomogeneous O - Z equation is obtained if the bulk density in E q . 1.3 is replaced w i t h the single particle density pfa), and is moved inside the integral. This means that />(r ) is another unknown, 3 and an additional relation between the single particle density and the pair correlation function is required to determine the variables completely. A l t h o u g h for liquids near a phase transition improved results are obtained from the inhomogenous equations, the implementation of the technique is significantly more difficult. B o t h versions of the O - Z equations are definitions for the direct correlation function, and are therefore formally exact. W i t h improved approximate closures, the homogeneous O - Z equation is expected to be useful lor systems near phase transitions, while maintaining computational simplicity [14]. In this thesis, the radial distribution functions of homogeneous liquids in contact with a hard macroparticle are calculated using the homogeneous O - Z integral equation with the H N C closure. Chapter 1. Introduction 6 As indicated previously, the thermodynamic quantities of the system can be evaluated from integrals over the radial distribution function. In particular, we calculate the excess energy, H e l m h o l t z free energy, and entropy of solvation of a single macroparticle. Here, the term "excess" refers to the difference between the thermodynamic quantity of a mixture and that of the pure solvent. In addition, we separate the excess quantities into two components; one due to the insertion of a Unite size particle into a l i q u i d , and the other due to the rearrangement of a solvent near a hard surface. T h e latter contribution is then associated with surface excess thermodynamic quantities. To obtain the excess internal energy at infinite dilution, we cannot simply evaluate the sum of pair energies over all solute-solvent interactions. T h e presence of a single solute particle perturbs the solvent. T h e structural changes in the solvent give rise to an additional contribution to the internal energy, which can be expressed as an integral over the derivative with respect to solute density of the solvent-solvent radial distribution function. Following Garisto ct.al. [15], the O - Z equation for the derivatives of the distribution functions are derived, and solved with the H N C closure. Thereby, the excess internal energy at infinite dilution is obtained. T h e excess Helmholtz free energy associated with adding a single particle to the bulk fluid is evaluated using a coupling integration technique [7]. A coupling parameter is defined in such a way that by varying it the solute is effectively added to or removed from solution. For a hard solute, the infinite dilution limit of the free energy is evaluated directly by applying some of the ideas of scaled particle theory [17]. T h e coupling parameter which multiplies the diameter of the hard particle, physically scales in a single solute. T h e change in free energy due to this scaling process is the quantity sought. A n expression for the excess Helmholtz free energy of a mixture is obtained, and subsequently, following Karplus et.al. [18], its infinite dilution limit is evaluated. Finally, the excess entropy of solvation is obtained from the thermodynamic expression relating the Chapter 1. Introduction 7 internal energy, Helmholtz free energy, and entropy. T h e remainder of this thesis is organized as follows. In Chapter 2, a detailed derivation of the equations required to calculate the density profiles and thermodynamic quantities are presented. Chapter 3, contains results obtained for a hard sphere solute in a hard sphere or Lennard-.Iones solvent at a number of different solvent densities and temperatures. In Chapter 4, results and conclusions are summarized. Chapter 2 Theory 2.1 T h e R a d i a l D i s t r i b u t i o n Function We consider a system consisting of particles ( 1 , N ) with spatial coordinates ( r i , ...r/v)- A distribution function expresses the probability of various configurations of a subset (1, ...,n) of these particles. T h e thermodynamic relationships in this chapter are derived in the canonical, ( N V T ) ensemble. T h e choice of ensemble is somewhat arbitrary, as for the calculation of thermodynamic quantities the various ensembles are indistinguishable [6]. Hence, the ensemble is chosen on the basis of computational convenience of thermodynamic quantities. W e consider fluids of particles which interact through pairwise additive potentials. G i v e n an N-particle system, pairwise additivity is equivalent to the condition U (r ...,r ) N where UN(TI , r ^ ) u N (2.4) = 5Z«(r j) , t is the total potential energy, u(vij) is the pair potential, and the summation is over all pairs of particles. For such a system the probability that particle 1 is in volume element r / r i , particle 2 in dr? , and particle n in volume element dr , n irrespective of the positions of the remaining particles is given by [2] N\ f... f (N-n)l v exp[-[1U ]dr ...dr N Z n+l N (2.5) N where the configurational integral is defined as (2.6) Chapter 2. Theory 9 T h e factor /? equals j^f, where is the Boltzniann constant. It can be shown, that in a bulk homogeneous fluid, density, p. For an isotropic fluid, p^ is equivalent to the bulk depends only on the separation between particles 1 and 2, and will be denoted ^ ' ( r 2 ) , or when no confusion is possible, simply p( \r). 2 If 2 r the probability of every particle i being in dvi at r,, were independent of the probability of every other particle j being in dtj at r-j, then p^ would be the product of single particle distributions, p^K W h e n this independence condition does not hold, a correlation function can be defined as / ( , l ) > (n,-,r ) = /) ( 1 ) w (r )..y '(r )/ (ri,...,r ) . i (2.7) , ) 1 n n T h e function c / ' ( r i , r , J is a factor that measures how much p ^ ( r - i , r ) n n deviates from the product of independent probabilities. Taking n = 2 in E q . 2.7, and noting that in a fluid, g^(r) is again a function only of the interparticle separations 7- , we 12 obtain the radial distribution function < / ' ( r i ) , which lias a simple physical meaning; it 2 2 is the factor which multiplies the bulk density, to give the local density. T h a t is, given a particle at the origin the number of particles in a spherical shell of thickness dr and radius ?• centered at the origin is pg 2.2 Choice of M o d e l T h e dominant features of pair potentials of simple liquids are a short-range strong repulsion, and a long-range weak attraction. T h e repulsion which is an expression of the fact that two particles cannot completely overlap, is a dominant contribution to the shortrange order exhibited by dense fluids. T h e attractive forces increase the stability of the fluid. From quantum mehanical perturbation theory, we know that at large separations the rare gas pair potential has the form [3] . u(r) = -c /r 6 G -c /r s 8 + ... . (2.8) Chapter 2. Theory 10 T h e exact form of the repulsive forces is not known. A l t h o u g h , in principle, an exact pair potential can be derived using quantum mechanics, the problem is too difficult numerically to be practical. Instead, simple analytical expressions w i t h parameters which can be adjusted to fit experimental data are employed. T h e simplest of these is the hard sphere potential, (2.9) where d is the hard sphere diameter. The hard sphere potential does not include attractive forces. Since the B o l t z m a n n factor e.x\)(—fiUN{T )) N is either zero or infinity for the potential defined by E q . 2.9, temperature has no effect on the average configurational properties of a hard sphere fluid. Despite its simplicity, a hard sphere l i q u i d exhibits many of the structural properties of real liquids. T h e Lennard-Jones(L-J) potential, incorporates the asymptotic behaviour of the true potential, and replaces the repulsion with an inverse power law, to obtain the mathematically simple expression, (2.10) w i t h cr and e as the adjustable parameters. T h i s potential has been used extensively for the study of argon-like liquids. T h e underlying assumption that the total potential energy can be expressed as a sum of pairwise additive potentials means that all manybody interactions are ignored. For rare-gas atoms, the net effect of the dominant manybody forces is repulsive [2]. We also know from experiment that the true potential has a deeper well and a weaker tail than the L-.l function [2]. Nevertheless, probably as a result of appropriate cancellation of errors, the L-.J potential has been reasonably successful in accounting for man}' of the properties of argon-like liquids [2]. Chapter 2. 2.3 Theory 11 Thermodynamics In this section, thermodynamic quantities in terms of radial distribution functions for a two component system will be presented. For classical liquids, all thermodynamic quantities can be separated into ideal and non-ideal parts. T h e non-ideal term contains the contribution of interparticle interactions. T h e terms energy, H e l m h o l t z free energy, and entropy will refer only to this second term, as it is the only contribution that will be considered here. The term "excess" will be used in the context of differences from the pure solvent. In particular, we calculate the change in energy, Helmholtz free energy, and entropy associated with taking a single solute particle from the ideal gas phase into solution in a pure solvent [15]. For the derivation of the excess potential energy, we appeal to the physical interpretation of g(r). Consider a two component mixture consisting of a solvent v and a solute fi. A r o u n d a central molecule of type //, the number of particles of the same component in the volume element A-Kr dr is given by 2 Pv<jvv(r)Anr dr, 2 where p u = p\ y i p is the bulk density, and \ „ is the mole fraction, jf. (2.11) Therefore, the potenial energy due to all solvent-solvent interactions in the m i x t u r e is (2.12) The factor N„ accounts for all solvent particles and the | ensures that each interaction is counted only once. Similarly, the contributions to the total potential energy from solvent-solute and solute-solute interactions are (2.13) Chapter 2. Theory respectively. T h e sum of these contributions divided by N is the potential energy per particle for the mixture T h e first term in E q . 2.18 represents the sum of all solute-solvent interactions. T h e second term is the contribution to the excess energy from the change in solvent structure induced by the presence of the solute. T h e next quantity sought is the change in Helmholtz free energy associated with adding a single solute to a pure solvent. A n expression for the H e l m h o l t z free energy in terms of the radial distribution function can be formulated v i a a coupling parameter integration technique [7]. A particle can be formally added to the system if the energy is expressed in terms of a coupling parameter £, as in Chapter 2. Theory 13 T h e parameter £ is defined so that u(rij;0) = 0 and u(rij;a) = u(vij), where a is a number signifying complete coupling. T h e value of a depends on the way we choose to introduce the effect of the particle. A solute can be effectively added to or removed from solution by physically scaling in its diameter, or by increasing the magnitude of its interaction potential from zero. T h e configurâtional integral, i n terms of the coupling parameter is Z {0 N = J exp[-(W (Ç)]dr ...dr N 1 . N (2.20) From the fundamental relation connecting thermodynamics to canonical ensemble averages, we have the configurational Helmholtz free energy as a function of £ [7], 0<A(Q>=-\n(?0) . (2.21) In order to calculate the free energy in excess of the pure solvent, the derivative of A with respect to £ is integrated, according to 0<A">= [ Jo l d t d l i < m o( > . (2.22) A t this point, it is convenient to introduce the function y ^{i"i 0 = exp[/3u(r; £)]flW( ' r u which is continuous for a l l potentials. If the configurational integral is differentiated explicitly with respect to £, an expression for the excess free energy in terms of £ and g(r;£) is obtained, f j < T h e factor \/N lt ^ > = ~ J*W*P J~ ^M-0u*Ar\mv**(rUydr. (2.23) is i m p l i c i t l y included, since by coupling a single particle, the excess free energy per solute particle is computed. For hard solutes it is useful to define a coupling parameter which physically scales in a particle by increasing its diameter. T h i s allows us to easily add particles of various sizes, and examine the free energy as a function of increasing solute diameter. Following the Chapter 2. 14 Theory approach of Reiss et.al. [17], the coupling parameter is defined in terms of the interaction of the hard sphere solute with the solvent, (2.24) Here, d is the solvent diameter and d is the variable representing the increasing diameter v of the solute. T h e particle is completely coupled when d = d^ and £ - - (d u + d^)/2 . T h e definition above implies that when the solute is completely uncoupled, its diameter is —d . v W h i l e seemingly unphysical, this in fact is necessary, because the effect of the particle does not disappear when the diameter is zero. A n impenetrable point particle would still be present when d = 0, which would exclude solvent particles from a volume element equal to the volume of a solvent particle. W h e n d = —d„, the effect of the solute vanishes, and the particle is completely removed. Once the potential is defined explicitly as a function of the coupling parameter, the derivative in the integrand of E q . 2.23 can be evaluated. For the potential and coupling parameter defined in E q . 2.24, the expression to be differentiated is (2.25) T h e function in E q . 2.25 is equivalent to 1 — H(Ç — r ) , where / / ( £ — is the Heav- iside step function, and its derivative with respect to £ is the Dirac delta function, —8(£ — r). Therefore the integration with respect to 7- becomes t r i v i a l , and the one- dimensional integral (2.26) is obtained. To evaluate the free energy for a hard sphere solute, then, requires an integration over the contact values of </(/•), for solute particles of increasing size. Chapter 2.4 2. Theory Integral equations The final piece of information required for the calculation of thermodynamic quantities is the radial distribution function. In order to find an approximate g(r) Ornstein-Zernike (O-Z) integral equation with the hypernetted-chain 2.4.1 we solve the ( H N C ) closure. T h e Ornstein-Zernike equation In 1914, Ornstein and Zernike [2] introduced the direct correlation function, c(r), v i a the relation, T h e function /),(ri,r ) represents the total pair correlation between particles 1 and 2 , 2 and is known as the total pair correlation function. We can see that / i ( r i , r ) does indeed 2 contain the total pair correlation between two particles, since when g(ri,r-i) = 1, and particles 1 and 2 are uncorrelated, /*(ri,r-2) is zero. The physical meaning of E q . 2.28 becomes clear if it is expanded recursively to give Chapter 2. Theory 16 T h e total correlation function consists of direct and indirect parts. T h e direct correlation function represents the direct influence of particles 1 and 2 on each other. T h e integral terms in E q . 2.30 represent the influence particle 1 on 2 as propagated by an increasing number of intermediate particles. T h e O-Z equation can be viewed as the definition of c ( r i , r ) , which is exact for a 2 fluid with pairwise additive potentials. For an isotropic homogeneous fluid, the O-Z equation has a simpler form which, when generalized for a mixture becomes for the present calculations of thermodynamic quantities. In order to evaluate the internal energy at infinite dilution according to E q . 2.18. we need to find the derivative with respect to solute density of the solvent-solvent radial distribution function. T h i s derivative can be determined, if following Garisto et.al [15], the pair correlation functions in the O - Z equation are expanded in powers of solute density as, Chapter 2. Theory Substituting Eqs. 2.34-2.37 for h(r) and c{r) into E q . 2.28, and collecting terms of equal power in solute density gives the O-Z equations for the pair correlation functions at infinite dilution, and their derivatives with respect to p^. T h e equations are U p o n Fourier transforming Eqs. 2.38-2.40, the integral equations in r-space become the following linear equations in k-space [2], Chapter 2. Theory 18 anc <"/^ = P» : 1 , ] [0 oc {k) + ^ j — uv 1 PI/CM , (2.46) (Vi/^A,'J respectively. In the formulation of the above equations the fact that /i„ = h^ M is used. u Since Eqs. 2.44 - 2.46 tire derived in the grand canonical ensemble [2], they are valid under constant volume conditions. under constant pressure. Experiments however, are likely to be performed In order to approximate a system under constant pressure, Eqs. 2.44 - 2.46 are rederived with the added restriction that upon addition of solute, an amount of solvent is removed to keep the total density constant. T h a t is, p u = p — p^ where p is the constant total reduced density. T h e O-Z equations for a system at constant total density are, 1- pcii(k) and 6 " M ~ T^ÏÏW) T T - m + i-pcK(t) ' ( 2 ' 4 8 ) T h e O - Z equation for the correlation functions at infinite dilution are the same for the constant density and constant volume formulations. T h e equation determining the derivatives of the solvent-solvent correlation functions however, contains an extra term. Therefore, as expected, thermodynamic quantities for the two cases will be different. 2.4.2 Choice of Closure Since, the O-Z equations are defining identities and both the total and the direct pair correlation functions are unknown, closures relating li(r) to C(Ï'), and Shfr) to 6c(r) are required to determine these functions completely. Several such closures have been proposed. T h e Perçus-Yevick ( P - Y ) approximation, //.(/•) - c ( r ) = e x p [ - i « ( r ) ] ( / ( r ) - l (2.49) Chapter 2. Theory 19 was derived by P e r ç u s and Yeviek in 1958 [7], using effective field techniques. T h e same equation was derived by Stell in 1963 [7], through graphical analysis of g(r) and c(r). T h e P - Y approximation is important because it has an analytical solution for hard sphere fluids. T h e mean spherical approximation ( M S A ) , (r) = 0, g c(r) = -pu(r), r < d , r > d , (2.50) is superior to the P - Y equation for square well fluids. Furthermore, it yields analytical solutions for some systems such as simple electrolyte and polar liquid models [2]. For hard sphere fluids the P - Y and M S A approximations are equivalent. A l l distribution functions defined here can be expressed as many-dimensional integrals over the particle coordinates. E q . 2.5 is an example of such an expression for the n-particle density. These integrals can be represented by graphs, and manipulated according to the topological properties of the diagrams [2]. Through such graphical analysis, diagrammatic expansions of c(r), and h(r) have been derived. From these expansions, the exact equation, hij(r) = djir) - H (r)] aij - flu^r) + B (r) tJ (2.51) resulted [2]. Together with the O-Z equation, E q . 2.51 yields two equations in the three unknowns, c(r), h(r) and B(r). The last term represents a class of diagrams known as elementary or bridge diagrams. Through further graphical analysis a t h i r d relation between g{r) and B(r) can be derived. T h i s however, is of no practical use since it introduces an infinite series of complex diagrams. Therefore, in practice B(r) is approximated. B y setting B(r) = 0, the hypernetted chain ( H N C ) approximation is obtained. A generalization of the H N C equation, the reference hypernetted chain ( R H N C ) equation, approximates B(r) by the bridge functions of a reference system w i t h a short range Chapter 2. Theory 20 potential. For the present calculations, the H N C approximation is used. In order to obtain a closure relation for Eqs. 2.38-2.40, the distribution functions in the H N C equation are expanded in powers of density according to Eqs. 2.34-2.37. Collecting terms of equal power in p gives the closure relations tl One of the advantages of choosing the H N C approximation is that an analytical expression can be obtained for the Helmholtz free energy at infinite dilution [18]. T h e formulation begins w i t h the expression for the Helmholtz free energy for a mixture at finite concentration, Chapter 2. Theory 21 T h e derivation of E q . 2.54 is somewhat lengthy, and is therefore given in the A p p e n d i x . To evaluate the excess Helmholtz free energy of solvation per solute particle, we first write the difference in free energy per solute particle of the mixture and the pure solvent, To evaluate the infinite dilution l i m i t , h.(r) and c(r) are again expanded in powers of solute density. In the limit of p^ = 0, this gives Chapter 2. Theory 22 solute particle at infinite dilution is /j < A ex > (2.61) = Pu T h e significance of this expression comes from the fact that in its derivation no assumptions about the form of the solute-solvent pair potential was required. Hence, unlike E q . 2.26, which is applicable only to hard solutes, E q . 2.61 is valid for a general potential. 2.4.3 N u m e r i c a l Details T h e H N C integral equations are solved by a standard iterative technique known as the P i c a r d method [2]. A n initial guess is made for the values of c;j(r). For our calculations, this initial guess for the homogeneous fluid is the analytical P - Y solution of a hard sphere fluid at the appropriate density. Calculations for mixtures of increasing solute size have the direct correlation functions of the mixture with the nearest smaller solute as input. These correlation functions are fast Fourier transformed, and used to calculate from equations (2.38) and (2.39).The function i]ff(r) fj^\k) is recovered by calculating the inverse Fourier transform of iiff{k). T h e ?/{^(r) are substituted in the closure relation to obtain a new value for cfj. T h i s procedure is repeated until the values of cfj converge. In order to force convergence, each guess is a combination of the results of the two previous iterations, computed according to output J-i + ac\ 'output ' (2.62) A value of 0.90 is used for a. A similar iteration procedure is used for the solution of Eqs. 2.46 and 2.53. T h e input consists of the infinite dilution results and the derivatives 6cjj(r) for the solute with the nearest smaller diameter. During each iteration SJ9J points are calculated, with grid size of 0.01. Chapter 3 Results and Discussion 3.1 Input Parameters The two systems studied were hard sphere and Lennard-Jones solvents, each w i t h a single hard sphere solute particle. T h e parameters for which calculations were performed are given in Table 3.1 , where p* = pa* is the reduced density and T* = ksT/e is the reduced temperature. System P* HS/HS 0.80 HS/HS 0.60 HS/LJ 0.80 1.35 HS/LJ 0.60 1.35 HS/LJ 0.50 1.35 HS/LJ 0.80 1.10 Table 3.1: Input Parameters; H S / H S refers to a hard sphere solvent with a hard sphere solute; H S / L J refers to a Lennard-Jones solvent with a hard sphere solute. 23 Chapter 3.2 3. Results and Discussion 24 Radial Distribution Functions T h e functions relevant to the calculation of thermodynamic quantities are (/[^(r), and <5<7„„(r). To illustrate the major trends observed, a representative sample of the calculated distribution functions are plotted. T h e variable ?• in the figures is in units of solvent diameter, d. In F i g . 3.1 the solute-solvent radial distribution functions for solutions consisting of a single hard sphere solute in a hard sphere solvent are plotted. represent T h e three curves calculated for three different solute sizes. T h e solid, dotted, and dashed curves correspond to solutes with diameters one, ten, and one-hundred times the diameter of the solvent respectively. T h e clotted and dashed curves have very similar shapes, illustrating the fact that <y|, J(r) becomes independent of the diameter of the macroparticle 0 for sufficiently large solutes. T h e discontinuity in each distribution function occurs at solute-solvent contact, given by (c/„ + d^)/2\ however, for comparison purposes, the plots oi gl°](r) are shifted toward the origin along the r-axis by a distance equal to difference between solute-solvent and solvent-solvent contact. F i g . 3.2 shows the solvent-solvent radial distribution functions for solutions consisting of a single hard sphere solute in a L-.J liquid. T h e meaning of the curves is analogous to those in F i g . 3.1. A g a i n , it is clear that as the diameter of the macroparticle increases <7^(r) becomes independent of solute size. Once the diameter of the solute is fifteen times that of the solvent, the shape of ^{^(r) remains virtually unchanged as the diameter of the macroparticle is increased. T h e same trend was observed by Perram and W h i t e [8] in their calculations for hard sphere fluids, using the P - Y approximation. T h e fact that becomes independent of solute size indicates that a large particle is a good representation of a flat wall. We emphasize that this is true only for separations which are small compared to the solute diameter. Chapter 3. Results and Discussion 25 Figure 3.1: Solute-solvent radial distribution functions for H S / H S solutions; p* = 0.8. T h e solid, dotted and dashed curves correspond to solvent-solute diameter ratios of 1:1, 1:10, and 1:100, respectively. Chapter 3. Results and Discussion 26 Figure 3.2: Solute-solvent radial distribution functions for H S / L J solutions; p* = 0.6, T* = 1.35. T h e solid, dotted, and dashed curves correspond to solvent-solute diameter ratios of 1:1, 1:10, and 1:100, respectively. Chapter 3. Results and Discussion 27 If the distance between the solute and the solvent particles is of the order of magnitude of the solute diameter, finite size effects become significant. For completeness, the solvent-solvent distribution function for the pure L - J solvent is given in F i g . 3.3. T h e function </f°J(r) for the pure hard sphere solvent is, of course, identical to ^[^(r) for the hard sphere m i x t u r e w i t h equal solute and solvent diameters, given in F i g . 3.1 by the solid curve. In F i g . 3.4, derivatives of the solvent-solvent distribution functions are plotted for a hard sphere solvent containing hard sphere solutes of increasing size. T h e same derivatives for a L - J solvent are given in F i g . 3.5. Since 8g (r) uu scales as the volume of the solute particle, the function is divided by the cube of the diameter. A s for g]fy(r), the shape of 8g {r)l vl/ d* becomes independent of the magnitude of solute diameter as the macroparticle is increased. In order to extract the surface term which is expected to scale as the surface area d , the general form of the derivative is assumed to be 2 8g„„(r, d) = Sgi{r)d + 8g2{r)d J 2 T h a t is, we assume that the dependence of 8g (r,d) ui/ arated from that on particle separation r. 8g (r)/d 2 1/u . (3.63) on solute diameter d, can be sep- In order to test the validity of E q . 3.63, for a fixed r as a function of solute diameter is plotted in Figs. 3.6-3.9. In the range where the assumption is valid, these plots should be linear. It is seen from Figs. 3.6-3.9 that at short distances from the macroparticle, the graph becomes linear for solutes as small as two to three times the solvent. As expected, for larger r , linearity is observed only for larger solutes. Given the values of 8g (r,d) for two different solute diameters, two equations in two ul/ unknowns are obtained, which can be solved to yield 8g\{r) and 8gw{r,di) = 8g\(r)d{ 8g ii\d ) = 8g\{r)d\ + 8g2{r)d\. vv 2 Sg2(r), + 8g2{r)d\ (3.64) Figure 3.3: Solvent-solvent radial distribution function for H S / L J solution; p* = 0.6, T* = 1.35. Figure 3.4: Derivatives of solvent-solvent radial distribution functions for H S / H S solutions; p* = 0.8. T h e solid, dotted, and dashed curves correspond to solvent-solute diameter ratios of 1:1, 1:10, and 1:100, respectively. Chapter 3. Results and Discussion 30 Figure 3.5: Derivatives of solvent-solvent radial distribution functions for H S / L J solutions; p* = 0.6, T* — 1.35. The solid, dotted, and dashed curves correspond to solvent-solute diameter ratios of 1:1, 1:10, and 1:100, respectively. Chapter 3. Results and Figure 3.6: 6g {r) ui/ Discussion at r = id for the H S / L J solution; p* = 0.6, T* = 1.35. 31 Chapter 3. Results and Discussion Chapter 3. Results and Discussion Chapter 3. Results and Figure 3.9: 8g (r) vv Discussion at r = 307/ for the H S / L J solution; p* = 0.6, T* = 1.35. Chapter 3. Results and Discussion 35 T h e functions obtained by solving Eqs. 3.64 using three different sets of solute diameters are plotted in Figs. 3.10-3.11. Specifically, Eqs. 3.64 were solved for the following macroparticle sizes d\ = I40d i d-2 = I20d i d = I40d i t d-2 = I00d i d\ = V20d i t d = 90d [ so vent x ao ven ao ven 2 so vent so vent sû vent . Figs. 3.10 and 3.11 show 8g\(r) and Sg2(r) for the full range of ?• at which the functions were calculated. The three curves obtained by solving Eqs. 3.64 using three different sets of solute diameters are identical up to a solute-solvent separation of approximately 20 solvent diameters. For r greater than 20d, the curves are considerably different. Therefore, given the macroparticle diameters used to calculate Sgl(r) and Sg2(r), the separation of surface and volume terms according to E q . 3.63 is valid up to solute-solvent distances of approximately 20d. The range of validity can be increased by performing calculations for larger solutes, but for the present calculations 20ri is more than sufficient. In Figs. 3.12 and 3.13, Sgi(r) and Sg2(r) are plotted on an expanded scale to show the short-range oscillatory structure of the functions. R e t u r n i n g now to E q . 2.40, the mathematical origin of the volume term can be identified. T h e term involving the solute-solvent correlation functions necessarily introduces a volume dependence due to the reduction of configuration space. T h i s dependence on the volume of the macroparticle is the result of inserting a particle of finite volume into solution. B y identifying the Sgl(r), term 8g2(r), the finite size effect is accounted for and the surface which contains the surface effect is extracted. F i n a l l y in F i g 3.14 the derivatives of the solvent-solvent distribution functions calculated w i t h the added restriction that the total density is constant are given. Chapter 3. Results and Discussion 36 Figure 3.10: Sgl(r) for H S / L J solution; p' = 0.6, T* = 1.35. T h e solid, dotted, and dashed lines correspond to solutions of Eqs. 3.61 using different solute diameter pairs. Chapter 3. Results and Discussion 37 Figure 3.11: 8g2{r) for H S / L J solution; p" = 0.6, T* = 1.35. T h e solid, dotted, and dashed lines correspond to solutions of Eqs. 3.61 using different solute diameter pairs. Chapter 3. Results and Discussion Figure 3.12: 8g\[r) for H S / L J solution on an expanded scale; p* = 0.6,T* = 1. Chapter 3. Results and Discussion Figure 3.13: 6g2(r) for H S / L J solution on an expanded scale; p* = 0.6, T* = 1. Chapter 3. Results and Discussion 40 Figure 3.14: 8g„„(r) at constant total density for H S / L J solution; p = 0.6, T* = 1.35. The solid, dashed, and dotted curves correspond to solvent-solute diameter ratios of 1:1, 1:10, and 1:100, respectively. Chapter 3. Results and 41 Discussion A s expected, the functions are qualitatively different from those calculated at constant volume (cf. F i g 3.2), but the trend that 8g„ (r)/d i u diameter for large solutes is again observed. becomes independent of solute Chapter 3.3 3. Results and Discussion 42 Thermodynamics T h e excess energy, entropy, and Helmholtz free energy of solvation per particle were calculated for each system with increasing solute, diameter. Figs. 3.15-3.17 depict these thermodynamic, quantities as a function of the square of the solute diameter. T h e choice of the function to be plotted is again dictated by the aim to separate the terms which are cubic i n solute diameter from those which are quadratic. In other words, it is assumed that the thermodynamic functions are of the form, A > < U > N, < S > < ex ex ex + A d' sur = u <P + U = 8 (P + S d d sur uol vot sur 2 2 2 . (3.65) T h e slope and intercept i n Figs. 3.15-3.17 give the volume and surface terms explicitly. Physically, the terms scaling as the volume and surface of the solute are associated with finite size effects and surface excess quantities respectively. M a t h e m a t i c a l l y it is clear that these terms should separate for the internal energy. In the. previous section it was shown that at short range 8 g „,,{!') can be separated into two terms; one which scales as the volume and the the other as the surface area of the macroparticle. For the systems studied here, 8g (r) ul/ completely determines the value of the internal energy, as the first term i n the energy expression E q . 2.18, is zero for hard sphere solutes. Also, since the major contribution of the potential is at short range, the short range structure of 8g (r) l/u determines the form of the internal energy. If the solute did not interact with the solvent through a hard potential, </[^](r) would contribute to the energy of solvation according to E q . 2.18. A s shown earlier, <7^(r) becomes independent of the diameter of sufficiently large solutes. Chapter 3. Results and Discussion Figure 3.15: Excess Energy of Solvation as a Function of Solute Diameter. Figure 3.16: Excess Free Energy of Solvation as a Function of Solute Diameter. Chapter 3. Results and Discussion Figure 3.17: Excess Entropy of Solvation as a Function of Solute Diameter. 45 Chapter 3. Results and Discussion 46 T h e only effect of the solute is to shift the function along the r axis, so that it becomes non-zero for 7" greater than solute-solvent contact. Since g\°l(r) is multiplied by r and 2 integrated over ?•, the shift introduces a t h i r d order dependence on the position of contact, known generally as the excluded volume effect. A surface dependence would arise from the term i n g\°\{r), since the number of solute-solvent interactions is proportional to the surface area of a macroparticle. Tables 3.2 and 3.3 contain the values of the two contributions to the internal energy of the L-.J solution obtained v i a two routes. The internal energy of a hard sphere solution is identically zero. T h e numbers in Table 3.2 come from numerically integrating E q . 2.18 for a number of solutes and obtaining the slope and intercept by extrapolating the linear region of plots such as F i g . 3.15 for each system studied. In Table 3.3, the volume and surface terms are calculated by integrating the equations, U = U = vot sur J47rr dr{ {r)6gl(r)} 2 Ul/l/ J4irr dr{u ,(r)tg2(r)} 2 u . (3.66) T h e only true difference between the two methods of calculation is that when using Eqs. 3.66, graphical analysis is not needed. T h e discrepancy in the numbers obtained for the contributions to the excess internal energy by the two methods is a result of numerical inaccuracies. If the graph of < U ex > /A r / t is extrapolated from sufficiently high values of d, the slope and intercept approach the values obtained by integrating Eqs. 3.66. T h e volume term of the internal energy is negative for all systems studied, indicating that when configuration space is reduced, particles rearrange in a way that lowers energy. T h e magnitude of this component increases from a density of 0.50 to 0.60, but decreases substantially from a density of 0.60 to 0.80. A t high density, the repulsive part of the potential is dominant, and the fact that the volume term becomes less negative can be attributed largely to the increase in effective density when the solute is added. Chapter 3. Results and Discussion 47 Table 3.2: Excess Energy of Solvation from extrapolation of < U cx System />* T* p < U ex > /JV > /N^d vs. d plots. 2 P < U ex M > /N„ Volume Term Surface T e r m HS/HS 0.60 0.0000 0.0000 HS/HS 0.80 0.0000 0.0000 HS/LJ 0.50 1.35 -0.5446 -1.3334 HS/LJ 0.60 1.35 -0.7201 -0.5979 HS/LJ 0.80 1.35 -0.3563 2.2374 HS/LJ 0.80 1.10 -0.7236 2.2700 Table 3.3: Excess Energy of Solvation from integration of 6 g\(r)and System P* T* p < U ex > /N 8g2(v). p<U > ex fl Volume Term Surface T e r m HS/HS 0.60 0.0000 0.0000 HS/HS 0.80 0.0000 0.0000 HS/LJ 0.50 1.35 -0.5438 -1.3998 HS/LJ 0.60 1.35 -0.7189 -0.6938 HS/LJ 0.80 1.35 -0.3554 2.1816 HS/LJ 0.80 1.10 -0.7215 2.0572 Chapter S. Results and Discussion 48 T h i s conclusion is supported by the values obtained for the internal energy at constant total density, given in Tables 3.4-3.5. T h e magnitude of the volume term does not decrease as much at p = 0.80 in the latter case, because the effective density is not increased upon introduction of the solute. T h e slight decrease in the magnitude of the volume component of the internal energy of high density liquids calculated at constant total density implies that at high density the rearrangement of solvent particles is less effective i n lowering the energy. T h e surface excess energy increases with solvent density, becoming positive at p = 0.80. T h e same trend is observed in the numbers calculated at constant total density. T h e surface excess energy at high density is greater in the former case because the actual solvent density is higher. W h e n the bulk density is low, the increase of the number of particles near the surface results in a negative contribution to the internal energy. For a high density liquid, repulsive forces become dominant as the number of solvent particles near the surface is increased, resulting in a positive contribution to the excess internal energy of solvation T h e Helmholtz free energy is also calculated using two different equations. Both methods involve obtaining the free energy for a number of solute particles with increasing diameters, and extrapolating the linear region of plots such as F i g . 3.16 to find the slope and intercept. T h e values in Table 3.6 were calculated by integrating E q . 2.26, while those in Table 3.7 come from using E q . 2.61. T h e latter has the advantage that it is not restricted to hard solutes. Both the volume and the surface contributions to the Helmholtz free energy are found to be positive. For a given density, the contribution of the volume term is greater in a hard sphere than in an L - J solvent. T h e component which scales as the surface area is smaller for a hard sphere fluid at low density and comparable for the two solvents at high density. Chapter 3. Results and Discussion 49 Table 3.4: Excess Energy of Solvation at constant total density from extrapolation of < U > /N^d vs. d plots. ex 2 System P* T* p < U tx > IN» p < U ex > /N^ Volume T e r m Surface T e r m HS/LJ 0.50 1.35 -0.5437 -1.4158 HS/LJ 0.60 1.35 -0.7261 -0.6005 HS/LJ 0.80 1.35 -0.6292 0.7689 Table 3.5: Excess Energy of Solvation at constant total density, from integration of 6gl(r)and 8g2(v). System P* T* p < U ex > 1 p < U ex > /N„ Volume T e r m Surface Term HS/LJ 0.50 1.35 -0.5442 -1.3626 HS/LJ 0.60 1.35 -0.7191 -0.6546 HS/LJ 0.80 1.35 -0.6284 0.7104 Chapter 3. Results and Discussion 50 Table .3.6: Excess Free Energy of Solvation, from E q . 2.26. System />* T* /j < A ex > /.\„ P < A ex > /7V Volume T e r m Surface T e r m HS/HS 0.60 1.5783 2.2080 HS/HS 0.80 4.0954 5.3757 HS/LJ 0.50 1.35 0.2253 1.1695 HS/LJ 0.60 1.35 0.4845 2.8352 HS/LJ 0.80 1.35 2.06927 7.2559 HS/LJ 0.80 1.10 1.89931 4.2404 M Table• 3.7: Excess Free Energy of Solvation, from E q . 2.61. System P* T* < A ex > (3 < A tx > /N„ Volume T e r m Surface Term HS/HS 0.60 1.6483 2.9238 HS/HS 0.80 4.2935 5.2179 HS/LJ 0.50 1.35 0.2342 1.4730 HS/LJ 0.60 1.35 0.5071 2.6999 HS/LJ 0.80 1.35 2.1667 6.4010 HS/LJ 0.80 1.10 1.9824 6.6915 Chapter 3. Results and Discussion 51 T h e values obtained for the components of the Helmholtz free energy by the two methods are considerably different. The two routes to the free energy are not thermodynamically different, and within the H N C approximation are mathematically equivalent as well. We believe that the source of the discrepancy in the numbers obtained by the two methods is numerical. T h e difference in the Helmholtz free energy for each solute between the two techniques is approximately 4-5 %. T h e discrepancy is magnified when the numbers are used in graphical analysis to determine the slope and intercept. In Table 3.8, calculated values for the excess entropy or solvation are given. These are obtained by using the thermodynamic relation, given in E q . 2.27, and subsequently finding the slope and intercept of graphs such as F i g . 3.17. B o t h contributions to the entropy are found to be negative, and becoming more negative w i t h increasing solvent density. T h e fact that the volume term is negative indicates that when configuration space is reduced with or without a corresponding increase in bulk density, the entropy of the system decreases. T h e negative contribution of the surface excess entropy is interpreted as the increased spatial ordering of the liquid near the surface. Chapter 3. Results and Discussion 52 Table 3.8: Excess Entropy of Solvation. System ^* T* fi < .S'"' > /k N„ B 0 < S* > e /k N„ B Volume Term Surface Term HS/HS 0.60 -1.5783 -2.2080 HS/HS 0.80 -4.0954 -5.3757 HS/LJ 0.50 1.35 -0.7699 -2.5029 HS/LJ 0.60 1.35 -1.2046 -2.8912 HS/LJ 0.80 1.35 -2.4256 -5.0185 HS/LJ 0.80 1.10 -2.70438 -4.5741 Chapter 4 Conclusions T h e Ornstein-Zernike ( O - Z ) equation for a two component m i x t u r e consisting of particles w i t h spherically symmetric potentials was solved with the hypernetted-chain ( H N C ) approximation to yield solute-solvent and solvent-solvent pair correlation functions at infinite dilution. T h e O - Z equation for the derivative of the solvent-solvent radial distribution function êg (r), is obtained by expressing the pair correlation functions as power ul/ series in solvent density. T h e function 6 </„„(»•) contains information about the rearrangement of solvent particles upon introduction of a solute particle. Solutions to both sets of coupled equations were obtained v i a an iterative procedure. Results were presented for solutes with diameters ranging in size from 1 to 150 times that of the solvent. T h e mixtures studied consist of hard solutes of varying sizes i n hard sphere and Lennard-Jones ( L - J ) solvents, at constant volume and a number of temperatures and solvent densities. A set of results was also presented for mixtures at constant total density. Each set of conditions defines a state well away from phase transitions. It is well known [9] that for large solutes the shape of the solute-solvent radial distribution function (r), becomes independent of the macroparticle diameter. Therefore, the large solute limit of a mixture at infinite dilution can be used to represent a bulk fluid in contact with a flat surface. In this thesis, the dependence of Sg„„(r) on solute diameter was investigated. T h e function Sg (r) ul/ was shown to be a sum of two dominant terms; one cubic 53 Chapter 4. Conclusions 54 and one quadratic in solute diameter. T h e two terms were separated, both graphically and algebraically. T h e cubic term was associated with the change in solvent-solvent correlations due to the exclusion from configuration space of the volume of the solute particle. It was found to be a function which has a short-range oscillatory structure and which decays to zero only at long-range. The term which scales as the surface area is interpreted as the function describing solvent rearrangement upon insertion of a flat surface into the pure solvent. T h e latter function also exhibits short-range oscillations, and contains a negative contribution to the total change in solvent correlations. T h e radial distribution functions and their derivatives were used to evaluate excess thermodynamic quantities of solvation. T h e excess internal energy, Helmholtz free energy, and entropy of solvation were shown to lie separable into cubic and quadratic terms in solute diameter. T h e excess internal energy of solvation at infinite dilution can be obtained from integrals over the solvent-solute radial distribution functions and the derivatives of the solvent-solvent radial distribution functions. In the present calculations, only hard solutes were considered. Hence, the excess internal energy of solvation is determined entirely by the change in solvent structure due to the presence of the solute, as the contribution of solvent-solute correlations is zero. A t low densities, both the volume and the surface terms of the internal energy were found to be negative. A t higher density, the volume term becomes more positive. T h i s effect is less pronounced when the total density is kept constant, and was therefore attributed in part to the increase in effective density upon insertion of a macroparticle into a liquid at constant volume. A s the density is increased, the surface excess energy also becomes more positive, and is greater than zero at high density. T h i s trend was interpreted as a result of increased solvent density near the surface. T h e excess Helmholtz free energy of solvation was calculated v i a two routes. One Chapter 4. Conclusions 55 technique is exact for an exact distribution function, but is restricted to hard solutes. T h e second method, which relies on the H N C approximation, is applicable to general solutes. T h e excess free energy was found to be positive for all solutions considered, indicating that inserting a macrosphere into a hard sphere or a L - J fluid increases the chemical potential of the liquid. From < U ex > /7V and < A ex M > jN^, the excess entropy of solvation was found by evaluating the thermodynamic relation connecting internal energy, Helmholtz free energy, and entropy. Both the surface and volume terms of the calculated excess entropies are negative, indicating that both finite size and surface el feet s involve increased spatial ordering of the liquid. T h e present approach can be extended to charged systems and liquids w i t h orientationally dependent potentials, such as water. In particular, contribution of orientational ordering to the entropy of liquids or liquid crystals near a flat surface can be studied. In its current form, the procedure can be used to obtain information about liquids interacting w i t h macroparticles and surfaces through hard and soft potentials. Bibliography [1] J . L . L e b o w i t z , E . M . W a i s m a n , The Liquid State of Matter: Fluids Simple and plex E . W . M o n t r o l l and J . L . L e b o w i t z , éd., N o r t h Holland Publishing, (1982) Com- [2] J - P . Hansen, I. R. M c D o n a l d , Theory of Simple Liquids, 2nd. ed., A c a d e m i c Press (1986) [3] D . Henderson, M . Plischke, J . Phys. C h e m . 92, 7177 (1988) [4] R . Zwanzig, J . C h e m . Phys. 22, 1420 (1954) [5] L . E . Reichl, Press,(1980) A Modern [6] T . L . H i l l , Statistical Course Mechanics, [7] D . A . M c Q u a r r i e , Statistical in Statistical Physics, University of Texas New York: Dover Publications Inc., (1987) Mechanics, Harper and Row, (1976) [8] J . W . P e r r a m , L . R. W h i t e , Faraday Discussions o the Chemical Society, 59, 29 (1975) [9] D . Henderson, J . A . Barker, F . F . A b r a h a m , M o l . Phys. 31, 1291 (1976) [10] L . S. S m i t h , L . L . Lee, J . C h e m . Phys. 71, 4085,(1979) [11] G . M . Torrie, P. G . K u s a l i k , G . N . Patey, J . C h e m . Phys. 88, 7826, (1988) [12] D . R . Berard, G . N . Patey, J . C h e m . Phys. 95, 5281 (1991) [13] R . Kjellander, S. Sarman, Molec. Phys. 70, 215 (1990) [14] D . R . Berard, G . N . Patey, (to be published) [15] F . Garisto, P. G . K u s a l i k , G . N . Patey, J . C h e m . Phys. 79, 6294 (1983) [16] S. Singer, D . Chandler, M o l . Phys. 55, 621 [17] L . Frisch, J . L Lebowitz, H . Reins. J . C h e m . Phys. 31, 369 (1959) [18] H - A . Y u , M . K a r p l u s , J . C h e m . Phys. 89, 2366, (1988) [19] F . L a d o , M o l . Phys. 47,283 (1982) 56 Appendix H e l m h o l t z Free E n e r g y for a M i x t u r e In this section, the H e l m h o l t z free energy for a two component m i x t u r e at finite concentration is derived. T h e details closely follow the analysis of fluids with non-spherical potentials by Lado [19], and the work of Singer ct.al. [16] with the extended RISM integral equations. To evaluate the H e l m h o l t z free energy arising from the interaction of particles, the coupling parameter integration method is again followed. T h i s time, however, the parameter, £, turns on and off all interparticle interactions. The total potential of the system is defined as a function of £ by LI {ri,r ; N iV 0 = £ "iA'u) • (A.l) i<j T h e subscript is included in the pair potential to emphasize that in a mixture the type of interaction depends upon the nature of particles /' and j . To continue, the H e l m h o l t z free energy as defined in the canonical ensemble is evaluated. T h e configuration integral now has the form, Z (Ç) N = J exp[-!UrH(r ...,r ;Z))dr . N v u N (A.2) A s Eqs. 2.21 and 2.22 indicate, the derivative of the Helmholtz free energy w i t h respect to the coupling parameter has to be found. W i t h the above configuration integral, this becomes, dp < A{Ç) > 1 dZ (Q --zW)~lîT57 N - (A 3) Evaluating E q . A . 3 explicitly yields, which is the canonical average of the potential energy. T h e above equation, therefore, can be written in terms of the radial distribuiou function as, The notation of E q . A . 5 is specific to an isotropic, homogeneous, two component mixture. The integrand can now be manipulated to a form which is an exact differential with respect to £. To achieve this, it is convenient to introduce two new functions, 59 U p o n rearrangement of E q . A.10 the integrand of E q . A . 8 can be expressed as a sum of an exact differential, and another term which needs to be manipulated further. The expression in its useful form is df d / 1 . dc 7 2 (A.ll) T h e total contribution of the second term to the Helmholtz free energy per particle of a m i x t u r e is, 2^Ex.-.Yi \\h%Vàr. (A.12) This s u m can be expressed as an exact differential of £ if O-Z equation for a mixture in fe-space is written in a m a t r i x form, h = c + h •p •c , (A.13) where (A.14) c (A.15) = and Pu 0 0 (A.16) P» M u l t i p l y i n g both sides of E q . A.13 by the matrix p, and adding the identity matrix gives, after some rearrangement, the relation, [I + ph]-' = [ I - p c ] . (A.17) we now have the integrand of A . 8 as an exact deferential with respect to the coupling parameter. Therefore, the integration with respect to £ is t r i v i a l . T h e final form of the Helmholtz free energy of a two component mixture is
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Statistical mechanics of solvation of macroparticles
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Statistical mechanics of solvation of macroparticles Vasarhelyi , Reka Z. 1992
pdf
Page Metadata
Item Metadata
Title | Statistical mechanics of solvation of macroparticles |
Creator |
Vasarhelyi , Reka Z. |
Date Issued | 1992 |
Description | Integral equation theories are employed to study the solvation of large particles. The solutions studied consist of spherically symmetric solvent particles which interact through hard sphere and Lennard-Jones potentials, and hard sphere solutes. In particular, the Ornstein-Zernike (O-Z) equation is solved with the hypernetted-chain (HNC) closure, to obtain the pair correlation functions of mixtures at infinite dilution. The pair correlation functions in the O-Z equation and the H N C closure are expressed as power series in solute density to yield a pair of coupled equations which determine the derivatives with respect to solute density of the solvent-solvent pair correlation functions. The latter describe the perturbation of the solvent upon addition of a single solute particle. The derivatives are analysed to yield components that scale as the volume and surface area of the macroparticle, and which are then identified as changes in solvent structure due to the presence of a finite size particle and a flat surface respectively. From the pair correlation functions and their derivatives the excess internal energy, Helmholtz free energy, and entropy of solvation are calculated. The excess quantities are also separated into contributions from finite size and surface effects. Both components of the excess internal energy are negative at low densities, and become less negative for high density liquids. The magnitude and sign of the two contributions to the energy depend on physical conditions such as temperature and pressure. The excess entropy of solvation is negative for all systems studied, indicating that introduction of a macroparticle or a flat surface increases the spatial ordering of a bulk liquid. |
Extent | 1458389 bytes |
Genre |
Thesis/Dissertation |
Type |
Text |
File Format | application/pdf |
Language | eng |
Date Available | 2008-12-18 |
Provider | Vancouver : University of British Columbia Library |
Rights | For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use. |
DOI | 10.14288/1.0061739 |
URI | http://hdl.handle.net/2429/3131 |
Degree |
Master of Science - MSc |
Program |
Chemistry |
Affiliation |
Science, Faculty of Chemistry, Department of |
Degree Grantor | University of British Columbia |
Graduation Date | 1992-11 |
Campus |
UBCV |
Scholarly Level | Graduate |
Aggregated Source Repository | DSpace |
Download
- Media
- 831-ubc_1992_fall_vasarhelyi_reka_z.pdf [ 1.39MB ]
- Metadata
- JSON: 831-1.0061739.json
- JSON-LD: 831-1.0061739-ld.json
- RDF/XML (Pretty): 831-1.0061739-rdf.xml
- RDF/JSON: 831-1.0061739-rdf.json
- Turtle: 831-1.0061739-turtle.txt
- N-Triples: 831-1.0061739-rdf-ntriples.txt
- Original Record: 831-1.0061739-source.json
- Full Text
- 831-1.0061739-fulltext.txt
- Citation
- 831-1.0061739.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.831.1-0061739/manifest