Electronic excitations in a dielectric continuum solvent with quantum Monte Carlo: Acrolein in water Franca Maria Floris, Claudia Filippi, and Claudio Amovilli Citation: The Journal of Chemical Physics 140, 034109 (2014); doi: 10.1063/1.4861429 View online: http://dx.doi.org/10.1063/1.4861429 View Table of Contents: http://scitation.aip.org/content/aip/journal/jcp/140/3?ver=pdfcov Published by the AIP Publishing Articles you may be interested in Electronic and vibrational spectra of protonated benzaldehyde-water clusters, [BZ-(H2O)n≤5]H+: Evidence for ground-state proton transfer to solvent for n ≥ 3 J. Chem. Phys. 140, 124314 (2014); 10.1063/1.4869341 Solvatochromic shifts of polar and non-polar molecules in ambient and supercritical water: A sequential quantum mechanics/molecular mechanics study including solute-solvent electron exchange-correlation J. Chem. Phys. 137, 214504 (2012); 10.1063/1.4769124 Quantum Monte Carlo formulation of volume polarization in dielectric continuum theory J. Chem. Phys. 129, 244106 (2008); 10.1063/1.3043804 Excitations in photoactive molecules from quantum Monte Carlo J. Chem. Phys. 121, 5836 (2004); 10.1063/1.1777212 Solvent-induced electronic decoherence: Configuration dependent dissipative evolution for solvated electron systems J. Chem. Phys. 116, 8429 (2002); 10.1063/1.1468887
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.230.73.202 On: Sun, 21 Dec 2014 17:45:36
THE JOURNAL OF CHEMICAL PHYSICS 140, 034109 (2014)
Electronic excitations in a dielectric continuum solvent with quantum Monte Carlo: Acrolein in water Franca Maria Floris,1,a) Claudia Filippi,2 and Claudio Amovilli1 1
Dipartimento di Chimica e Chimica Industriale, Università di Pisa, Via Risorgimento 35, 56126 Pisa, Italy MESA+ Institute for Nanotechnology, University of Twente, P.O. Box 217, 7500 AE Enschede, The Netherlands 2
(Received 9 September 2013; accepted 23 December 2013; published online 17 January 2014) We investigate here the vertical n → π ∗ and π → π ∗ transitions of s-trans-acrolein in aqueous solution by means of a polarizable continuum model (PCM) we have developed for the treatment of the solute at the quantum Monte Carlo (QMC) level of the theory. We employ the QMC approach which allows us to work with highly correlated electronic wave functions for both the solute ground and excited states and, to study the vertical transitions in the solvent, adopt the commonly used scheme of considering fast and slow dielectric polarization. To perform calculations in a non-equilibrium solvation regime for the solute excited state, we add a correction to the global dielectric polarization charge density, obtained self consistently with the solute ground-state wave function by assuming a linear-response scheme. For the solvent polarization in the field of the solute in the ground state, we use the static dielectric constant while, for the electronic dielectric polarization, we employ the solvent refractive index evaluated at the same frequency of the photon absorbed by the solute for the transition. This choice is shown to be better than adopting the most commonly used value of refractive index measured in the region of visible radiation. Our QMC calculations show that, for standard cavities, the solvatochromic shifts obtained with the PCM are underestimated, even though of the correct sign, for both transitions of acrolein in water. Only by reducing the size of the cavity to values where more than one electron is escaped to the solvent region, we regain the experimental shift for the n → π ∗ case and also improve considerably the shift for the π → π ∗ transition. © 2014 AIP Publishing LLC. [http://dx.doi.org/10.1063/1.4861429] I. INTRODUCTION
The understanding of solvent effects is of fundamental importance for the study of a great variety of processes involving molecules dispersed in a condensed phase.1 In chemical solutions, solvent effects can also influence the absorption of UV and visible radiation with effects that can be measured by spectroscopic techniques and exploited in photochemistry.2 In this work, we deal with the study of how an environment, like the one generated by a solvent, changes the electronic energy of a given reference system in a fast process such as the absorption of a photon. It is well established that this energy change originates from intermolecular forces but the complexity of the liquid phase, which characterizes all solutions, renders difficult a simple rationalization. The most widely used theoretical models of solvation start from the often valid approximation that the solute, that is, the reference system, maintains its individuality. With this important assumption, the environment can be added to the model at various levels of sophistication. The most natural way is to add a given, generally large number of solvent molecules.3 However, the geometrical distribution of these molecules cannot be automatically assigned, owing to the weak intermolecular interactions, and it must be determined through computer simulations a) Author to whom correspondence should be addressed. Electronic mail:
[email protected]
0021-9606/2014/140(3)/034109/9/$30.00
in order to achieve a statistically significant description.4 Moreover, the treatment of the solvent molecules can range from the same quantum mechanical level used for the solute5 to a classical mechanical level if reliable potential energy functions are available to describe intermolecular interactions.6 The employed number of explicit solvent molecules depends strongly on the level of treatment and can be large only for the simplest classical mechanics model. Computer simulations can then be performed by standard Monte Carlo (MC)7 or molecular dynamics (MD)8 methods. At the opposite, a completely different approach, which avoids heavy computer simulations, is the so-called polarizable continuum model (PCM).9 In this model, the solvent is represented by a continuum dielectric medium to describe the solute molecules in a cavity of a given shape.10 Dielectric continuum theory also allows one to approximately include the effects of the electrostatic potential produced by the polarization of the solvent in the presence of a solute. These effects are responsible for a significant fraction of the free energy of solvation, namely, the so-called polarization free energy.11 Even though the polarization free energy is generally not the leading term in the decomposition of the free energy of solvation,13 it is the most important contribution due to the solvent in processes occurring in solution including absorption of photons. Nowadays, there are numerous variants of such PCM, which also include Pauli repulsion and dispersion interactions.12, 13 In the PCM, the solute is treated at
140, 034109-1
© 2014 AIP Publishing LLC
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.230.73.202 On: Sun, 21 Dec 2014 17:45:36
034109-2
Floris, Filippi, and Amovilli
J. Chem. Phys. 140, 034109 (2014)
a quantum mechanical level through the use of a modified Hamiltonian containing a term representing the solvent reaction field.9, 11 The complete solution of the related physical problem is achieved by coupling the Poisson’s equation to the quantum mechanical description of the solute. Recently, we have proposed a PCM designed for the quantum Monte Carlo (QMC) treatment of the solute.14, 15 In this work, we present the extension of our model to the study of fast electronic transitions due to photon absorption in solution. This non-trivial extension is necessary to construct the appropriate environment for the non-equilibrium solvation of a solute in an electronic excited state in the framework of a surface and volume polarization for electrostatics (SVPE)16 refined method. As a test example, we consider the vertical n → π ∗ and π → π ∗ transitions of s-trans-acrolein in aqueous solution. We choose this system because it has been widely studied, it shows a solvent effect in the UV absorption spectrum and, for the two aforementioned transitions, the effect of the solvent water is opposite.17–20 The outline of the paper is as follows. In Sec. II, we review our QMC implementation of the PCM and, in Subsection II D, we illustrate the extension of the method to treat electronic vertical excitations of solutes. In Sec. III, we describe the details of the calculations performed in this work and, in Sec. IV, we report the QMC calculations with a discussion of results. Finally, we give a brief summary and discuss possible prospective for further work. II. THEORY
In this section, we review the theory we have developed in recent works,14, 15 where we implemented a PCM version specifically designed for QMC calculations. We then illustrate the necessary modifications to treat fast vertical excitations according to the common standard adopted for PCM in this kind of applications (see, for example, Ref. 21).
B. The solvent reaction field
For the solute in solution, the Schrödinger equation is coupled to the Poisson’s equation,11, 28 divE = 4πρ (tot) ,
(3) (tot)
is the total charge where E is the total electric field and ρ density given by the nuclear and electronic charge distribution of the solute, and by the polarization charge distribution of the dielectric induced by the solute. The solute molecule is placed within a cavity created in the dielectric medium of dielectric constant . Because the dielectric constant inside the cavity has the vacuum value, there is a discontinuity at the cavity surface, which determines a surface polarization charge density. Moreover, when the solute electrons escape from the cavity, we have also a volume polarization charge density defined in the solvent region. These two polarization charge densities produce the so-called solvent reaction field, which contributes to the total electric field in Eq. (3) above. In the most widely used PCM,29 there is not a volume polarization charge density and, generally, the surface polarization charge density includes corrections to account implicitly for the solute wave function tails outside of the cavity. In contrast, the PCM used in this work is a specific implementation for QMC calculations developed within the general approach known as the surface and volume polarization for electrostatics method.16 In the following, we first recall the main equations of the method presented in our previous works14, 15 and, then, illustrate the extensions of our method to treat vertical electronic transitions. 1. Volume polarization
The so-called volume polarization charge density is given
A. The solute
by
The QMC many-body wave function used here to describe the solute is written in a spin-assigned Slater-Jastrow form22–24 as (r1 , r2 , . . .) = (r1 , r2 , . . .)J (r1 , r2 , . . .) , with determinantal component given by ↑ ↓ = DK DK dK ,
(1)
(2)
K ↑
correlation energy with a relatively small number of determinants and can be a good trial wave function for a fixed-node diffusion Monte Carlo (DMC)22, 27 calculation.
↓
where DK and DK are the Slater determinants constructed from the occupied orbitals of spin-up and spin-down electrons, respectively, and dK are the mixing coefficients. The Jastrow correlation factor J is the exponential of the sum of three fifth-order polynomials of the electron-nuclear (e-n), the electron-electron (e-e), and of pure 3-body mixed e-e and e-n distances.25 The parameters in the Jastrow factor, the determinantal expansion, and the orbitals are optimized at the variational Monte Carlo (VMC) level in energy minimization.26 This wave function is able to include a substantial part of the
(pol)
ρvol (r) =
1 − 1 ρe (r) ,
(4)
where ρ e (r), in atomic units, is the electronic density of the solute. According to Eq. (4), the volume polarization charge density is non-zero only where differs from 1, the vacuum value of the dielectric constant. Therefore, the electrostatic potential associated with this charge density is given by the following integral over the domain outside the cavity C: N (pol) ρvol (r ) θ (ri ) 1 (pol) dr | | , = 1− φvol (r) = | |r − r |r − ri | r ∈ C i=1 (5) where is the solute wave function and θ (r) is simply given by θ (r) = 0 (r ∈ C) ,
θ (r) = 1 (r ∈ C) .
(6)
We employ the variational Monte Carlo approach to compute the integral above by sampling a set of M configurations
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.230.73.202 On: Sun, 21 Dec 2014 17:45:36
034109-3
Floris, Filippi, and Amovilli
J. Chem. Phys. 140, 034109 (2014)
(k) 2 (r(k) 1 , . . . , rN ) from the square of the wave function || using the Metropolis Monte Carlo method. The values of the integrand computed with these configurations are averaged to give N M 1 1 θ (r(k) (pol) i ) φvol (r) ≈ 1− , (7) M k=1 i=1 |r − r(k) i |
where the estimate tends to the exact result as M → ∞. By comparing this expression with the analogous equation for a discretized set of nc volume polarization charges, (pol)
φvol (r) ≈
nc l=1
ql , |r − rl |
(8)
we see that our expression corresponds to nc ≤ N × M point charges with charge ql = (1 − 1/)/M and positions corresponding to the sampled one-electron coordinates outside the (k) cavity, rl = r(k) i if θ (ri ) = 1.
Once the electric field has been obtained after a selfconsistent procedure, the surface polarization charge density σ can be computed at the position r on the cavity border, according to the relation (1 − ) n+ · E− (r) , (9) 4π where n+ is the unity vector pointing outside the cavity from r and E− is the total electric field evaluated at the surface immediately inside of the cavity. The electric field results from the superposition of the fields generated by free and polarization charges, and can be written as σ (r) =
(10)
Thus, at the end of a VMC run in the solvent, the normal component of the averaged solute field is computed as well as the corresponding quantities related to the fields produced by polarization charges. The normal component of the field produced by volume polarization is easily computed from Eq. (8) while, for the same quantity related to the surface polarization charge density, we use the relation n+ · Esurf− (r) = −2π σ (r) n+ · (r − ra ) σ (ra ) da , + |r − ra |3
σ (rk ) ≈
qk , a
(12)
with
2. Surface polarization
E = Esolute + Esurf + Evol .
a seam between two spheres, a point charge is assigned to a portion of surface with area a given by the inverse of the number density p. At the seams where two point charges of different spheres are closer than a fixed length proportional to √ a, the two point charges are replaced by a single one in the middle, and all properties related to them are averaged with the same weights. The coordinates of the point charges are separately determined on each sphere before these are joined to model the cavity. To obtain a homogeneous distribution of points on each sphere, we resort to the Thomson problem that we solve by numerical simulation.30, 31 In the relevant quadrature formula,14, 15, 31 each of these point charges is weighted by a factor to account for the discrete summation in the surface integral in Eq. (11) above. Thus, the surface charge density σ at the position rk of the point charge qk is approximated as
qk =
Gkj ( , ) n+ · Evol (rj ) + Esolute (rj ) ,
where Gkj depends only on the shape of the cavity and on the solvent dielectric constant. Equation (13) establishes a relation between the Poisson’s equation and the solute Schrödinger equation from which the electron density (and subsequently the volume polarization charges) is derived, and is the basis for a self-consistent procedure which leads to a simultaneous solution of both equations.
C. The ground-state wave function
First, all determinantal, Jastrow, and orbital parameters of the solute wave function are simultaneously optimized within VMC by energy minimization in vacuo. Then, a starting set of volume and surface polarization charges are obtained through Eqs. (7), (8), and (13) and inserted in the following energy functional: F [] = |Helec | +
α
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.230.73.202 On: Sun, 21 Dec 2014 17:45:36
THE JOURNAL OF CHEMICAL PHYSICS 140, 034109 (2014)
Electronic excitations in a dielectric continuum solvent with quantum Monte Carlo: Acrolein in water Franca Maria Floris,1,a) Claudia Filippi,2 and Claudio Amovilli1 1
Dipartimento di Chimica e Chimica Industriale, Università di Pisa, Via Risorgimento 35, 56126 Pisa, Italy MESA+ Institute for Nanotechnology, University of Twente, P.O. Box 217, 7500 AE Enschede, The Netherlands 2
(Received 9 September 2013; accepted 23 December 2013; published online 17 January 2014) We investigate here the vertical n → π ∗ and π → π ∗ transitions of s-trans-acrolein in aqueous solution by means of a polarizable continuum model (PCM) we have developed for the treatment of the solute at the quantum Monte Carlo (QMC) level of the theory. We employ the QMC approach which allows us to work with highly correlated electronic wave functions for both the solute ground and excited states and, to study the vertical transitions in the solvent, adopt the commonly used scheme of considering fast and slow dielectric polarization. To perform calculations in a non-equilibrium solvation regime for the solute excited state, we add a correction to the global dielectric polarization charge density, obtained self consistently with the solute ground-state wave function by assuming a linear-response scheme. For the solvent polarization in the field of the solute in the ground state, we use the static dielectric constant while, for the electronic dielectric polarization, we employ the solvent refractive index evaluated at the same frequency of the photon absorbed by the solute for the transition. This choice is shown to be better than adopting the most commonly used value of refractive index measured in the region of visible radiation. Our QMC calculations show that, for standard cavities, the solvatochromic shifts obtained with the PCM are underestimated, even though of the correct sign, for both transitions of acrolein in water. Only by reducing the size of the cavity to values where more than one electron is escaped to the solvent region, we regain the experimental shift for the n → π ∗ case and also improve considerably the shift for the π → π ∗ transition. © 2014 AIP Publishing LLC. [http://dx.doi.org/10.1063/1.4861429] I. INTRODUCTION
The understanding of solvent effects is of fundamental importance for the study of a great variety of processes involving molecules dispersed in a condensed phase.1 In chemical solutions, solvent effects can also influence the absorption of UV and visible radiation with effects that can be measured by spectroscopic techniques and exploited in photochemistry.2 In this work, we deal with the study of how an environment, like the one generated by a solvent, changes the electronic energy of a given reference system in a fast process such as the absorption of a photon. It is well established that this energy change originates from intermolecular forces but the complexity of the liquid phase, which characterizes all solutions, renders difficult a simple rationalization. The most widely used theoretical models of solvation start from the often valid approximation that the solute, that is, the reference system, maintains its individuality. With this important assumption, the environment can be added to the model at various levels of sophistication. The most natural way is to add a given, generally large number of solvent molecules.3 However, the geometrical distribution of these molecules cannot be automatically assigned, owing to the weak intermolecular interactions, and it must be determined through computer simulations a) Author to whom correspondence should be addressed. Electronic mail:
[email protected]
0021-9606/2014/140(3)/034109/9/$30.00
in order to achieve a statistically significant description.4 Moreover, the treatment of the solvent molecules can range from the same quantum mechanical level used for the solute5 to a classical mechanical level if reliable potential energy functions are available to describe intermolecular interactions.6 The employed number of explicit solvent molecules depends strongly on the level of treatment and can be large only for the simplest classical mechanics model. Computer simulations can then be performed by standard Monte Carlo (MC)7 or molecular dynamics (MD)8 methods. At the opposite, a completely different approach, which avoids heavy computer simulations, is the so-called polarizable continuum model (PCM).9 In this model, the solvent is represented by a continuum dielectric medium to describe the solute molecules in a cavity of a given shape.10 Dielectric continuum theory also allows one to approximately include the effects of the electrostatic potential produced by the polarization of the solvent in the presence of a solute. These effects are responsible for a significant fraction of the free energy of solvation, namely, the so-called polarization free energy.11 Even though the polarization free energy is generally not the leading term in the decomposition of the free energy of solvation,13 it is the most important contribution due to the solvent in processes occurring in solution including absorption of photons. Nowadays, there are numerous variants of such PCM, which also include Pauli repulsion and dispersion interactions.12, 13 In the PCM, the solute is treated at
140, 034109-1
© 2014 AIP Publishing LLC
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.230.73.202 On: Sun, 21 Dec 2014 17:45:36
034109-2
Floris, Filippi, and Amovilli
J. Chem. Phys. 140, 034109 (2014)
a quantum mechanical level through the use of a modified Hamiltonian containing a term representing the solvent reaction field.9, 11 The complete solution of the related physical problem is achieved by coupling the Poisson’s equation to the quantum mechanical description of the solute. Recently, we have proposed a PCM designed for the quantum Monte Carlo (QMC) treatment of the solute.14, 15 In this work, we present the extension of our model to the study of fast electronic transitions due to photon absorption in solution. This non-trivial extension is necessary to construct the appropriate environment for the non-equilibrium solvation of a solute in an electronic excited state in the framework of a surface and volume polarization for electrostatics (SVPE)16 refined method. As a test example, we consider the vertical n → π ∗ and π → π ∗ transitions of s-trans-acrolein in aqueous solution. We choose this system because it has been widely studied, it shows a solvent effect in the UV absorption spectrum and, for the two aforementioned transitions, the effect of the solvent water is opposite.17–20 The outline of the paper is as follows. In Sec. II, we review our QMC implementation of the PCM and, in Subsection II D, we illustrate the extension of the method to treat electronic vertical excitations of solutes. In Sec. III, we describe the details of the calculations performed in this work and, in Sec. IV, we report the QMC calculations with a discussion of results. Finally, we give a brief summary and discuss possible prospective for further work. II. THEORY
In this section, we review the theory we have developed in recent works,14, 15 where we implemented a PCM version specifically designed for QMC calculations. We then illustrate the necessary modifications to treat fast vertical excitations according to the common standard adopted for PCM in this kind of applications (see, for example, Ref. 21).
B. The solvent reaction field
For the solute in solution, the Schrödinger equation is coupled to the Poisson’s equation,11, 28 divE = 4πρ (tot) ,
(3) (tot)
is the total charge where E is the total electric field and ρ density given by the nuclear and electronic charge distribution of the solute, and by the polarization charge distribution of the dielectric induced by the solute. The solute molecule is placed within a cavity created in the dielectric medium of dielectric constant . Because the dielectric constant inside the cavity has the vacuum value, there is a discontinuity at the cavity surface, which determines a surface polarization charge density. Moreover, when the solute electrons escape from the cavity, we have also a volume polarization charge density defined in the solvent region. These two polarization charge densities produce the so-called solvent reaction field, which contributes to the total electric field in Eq. (3) above. In the most widely used PCM,29 there is not a volume polarization charge density and, generally, the surface polarization charge density includes corrections to account implicitly for the solute wave function tails outside of the cavity. In contrast, the PCM used in this work is a specific implementation for QMC calculations developed within the general approach known as the surface and volume polarization for electrostatics method.16 In the following, we first recall the main equations of the method presented in our previous works14, 15 and, then, illustrate the extensions of our method to treat vertical electronic transitions. 1. Volume polarization
The so-called volume polarization charge density is given
A. The solute
by
The QMC many-body wave function used here to describe the solute is written in a spin-assigned Slater-Jastrow form22–24 as (r1 , r2 , . . .) = (r1 , r2 , . . .)J (r1 , r2 , . . .) , with determinantal component given by ↑ ↓ = DK DK dK ,
(1)
(2)
K ↑
correlation energy with a relatively small number of determinants and can be a good trial wave function for a fixed-node diffusion Monte Carlo (DMC)22, 27 calculation.
↓
where DK and DK are the Slater determinants constructed from the occupied orbitals of spin-up and spin-down electrons, respectively, and dK are the mixing coefficients. The Jastrow correlation factor J is the exponential of the sum of three fifth-order polynomials of the electron-nuclear (e-n), the electron-electron (e-e), and of pure 3-body mixed e-e and e-n distances.25 The parameters in the Jastrow factor, the determinantal expansion, and the orbitals are optimized at the variational Monte Carlo (VMC) level in energy minimization.26 This wave function is able to include a substantial part of the
(pol)
ρvol (r) =
1 − 1 ρe (r) ,
(4)
where ρ e (r), in atomic units, is the electronic density of the solute. According to Eq. (4), the volume polarization charge density is non-zero only where differs from 1, the vacuum value of the dielectric constant. Therefore, the electrostatic potential associated with this charge density is given by the following integral over the domain outside the cavity C: N (pol) ρvol (r ) θ (ri ) 1 (pol) dr | | , = 1− φvol (r) = | |r − r |r − ri | r ∈ C i=1 (5) where is the solute wave function and θ (r) is simply given by θ (r) = 0 (r ∈ C) ,
θ (r) = 1 (r ∈ C) .
(6)
We employ the variational Monte Carlo approach to compute the integral above by sampling a set of M configurations
This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 131.230.73.202 On: Sun, 21 Dec 2014 17:45:36
034109-3
Floris, Filippi, and Amovilli
J. Chem. Phys. 140, 034109 (2014)
(k) 2 (r(k) 1 , . . . , rN ) from the square of the wave function || using the Metropolis Monte Carlo method. The values of the integrand computed with these configurations are averaged to give N M 1 1 θ (r(k) (pol) i ) φvol (r) ≈ 1− , (7) M k=1 i=1 |r − r(k) i |
where the estimate tends to the exact result as M → ∞. By comparing this expression with the analogous equation for a discretized set of nc volume polarization charges, (pol)
φvol (r) ≈
nc l=1
ql , |r − rl |
(8)
we see that our expression corresponds to nc ≤ N × M point charges with charge ql = (1 − 1/)/M and positions corresponding to the sampled one-electron coordinates outside the (k) cavity, rl = r(k) i if θ (ri ) = 1.
Once the electric field has been obtained after a selfconsistent procedure, the surface polarization charge density σ can be computed at the position r on the cavity border, according to the relation (1 − ) n+ · E− (r) , (9) 4π where n+ is the unity vector pointing outside the cavity from r and E− is the total electric field evaluated at the surface immediately inside of the cavity. The electric field results from the superposition of the fields generated by free and polarization charges, and can be written as σ (r) =
(10)
Thus, at the end of a VMC run in the solvent, the normal component of the averaged solute field is computed as well as the corresponding quantities related to the fields produced by polarization charges. The normal component of the field produced by volume polarization is easily computed from Eq. (8) while, for the same quantity related to the surface polarization charge density, we use the relation n+ · Esurf− (r) = −2π σ (r) n+ · (r − ra ) σ (ra ) da , + |r − ra |3
σ (rk ) ≈
qk , a
(12)
with
2. Surface polarization
E = Esolute + Esurf + Evol .
a seam between two spheres, a point charge is assigned to a portion of surface with area a given by the inverse of the number density p. At the seams where two point charges of different spheres are closer than a fixed length proportional to √ a, the two point charges are replaced by a single one in the middle, and all properties related to them are averaged with the same weights. The coordinates of the point charges are separately determined on each sphere before these are joined to model the cavity. To obtain a homogeneous distribution of points on each sphere, we resort to the Thomson problem that we solve by numerical simulation.30, 31 In the relevant quadrature formula,14, 15, 31 each of these point charges is weighted by a factor to account for the discrete summation in the surface integral in Eq. (11) above. Thus, the surface charge density σ at the position rk of the point charge qk is approximated as
qk =
Gkj ( , ) n+ · Evol (rj ) + Esolute (rj ) ,
where Gkj depends only on the shape of the cavity and on the solvent dielectric constant. Equation (13) establishes a relation between the Poisson’s equation and the solute Schrödinger equation from which the electron density (and subsequently the volume polarization charges) is derived, and is the basis for a self-consistent procedure which leads to a simultaneous solution of both equations.
C. The ground-state wave function
First, all determinantal, Jastrow, and orbital parameters of the solute wave function are simultaneously optimized within VMC by energy minimization in vacuo. Then, a starting set of volume and surface polarization charges are obtained through Eqs. (7), (8), and (13) and inserted in the following energy functional: F [] = |Helec | +
α
