Development of Reaction Density Functional Theory and Its

Subscriber access provided by University of South Dakota

C: Surfaces, Interfaces, Porous Materials, and Catalysis

Development of Reaction Density Functional Theory and Its Application to Glycine Tautomerization Reaction in Aqueous Solution Weiqiang Tang, Cheng Cai, Shuangliang Zhao, and Honglai Liu J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.8b05383 • Publication Date (Web): 22 Aug 2018 Downloaded from http://pubs.acs.org on August 23, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

1

Development of Reaction Density Functional Theory and Its Application to

2

Glycine Tautomerization Reaction in Aqueous Solution

3

Weiqiang Tang1, Cheng Cai2, Shuangliang Zhao1,1 and Honglai Liu2 1

4

State Key Laboratory of Chemical Engineering and School of Chemical Engineering, East

5

China University of Science and Technology, Shanghai, 200237, China 2

6

School of Chemistry and Molecular Engineering, East China University of Science and

7

Technology, Shanghai, 200237, China

8

Abstract

9

Accurate determination of reaction barrier is crucial in chemical reaction engineering since

10

it’s directly associated with the reaction rate and reaction path. Whereas massive quantum

11

calculations relying on commercial software have been reported for understanding reaction

12

mechanisms, the solvent effect on the reaction barrier is still poorly understood. The difficulty

13

originates partially from the ill capture of inhomogeneous solvent structure around the

14

reaction site. Herein, we propose a reaction density functional theory (RxDFT) by combining

15

quantum density functional theory for calculating intrinsic reaction energy with classical

16

density functional theory for addressing solvation contribution. For demonstration, the

17

glycine tautomerization reaction in aqueous solution is revisited by using RxDFT. We

18

illustrate that the free energy profiles of two distinct reaction paths can be well predicted, in

19

comparison with the results of MD/MP2 simulation. Since no solvent molecule is explicitly

20

considered, the developed RxDFT is promising to provide an efficient alternative to QM/MM

21

approach or ab initio calculation.

22

1

To whom correspondence should be addressed. Email: [email protected]

1

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

I.

Introduction

2

It has been well recognized that reaction solvent has a profound effect on the reaction rate,

3

chemical equilibrium and even the reaction mechanism.1-3 The use of solvents in organic

4

chemical reactions, especially in synthesis, has become a regulation method. Using different

5

solvents, the reaction efficiencies are quite different.4 For example, in the substitution reaction

6

of butyl bromide and sodium cyanide: C 4 H 9 Br + NaCN → C 4 H 9 CN + NaBr , the solvent plays

7

a pivotal role.5 Namely, if butyl bromide is simply mixed with sodium cyanide in aqueous

8

solution, the reaction does not occur even if refluxed at 373 K for two weeks. This is mainly

9

due to the fact that butyl bromide is insoluble in water and the substrate cannot contact with the

10

reagents. In contrast, the reaction can take place in the presence of alcoholic solvents although

11

the reaction rate and yield are still quite low. Interestingly, the reaction rate in

12

dimethylformamide (DMF) solution is ten times faster than that in alcohol. The influences of

13

the solvent on the reaction process are collectively referred to as the solvent effect or the

14

medium effect of the reaction.4

15

Generally, there are three categories of methods for addressing the solvent effects of

16

chemical reactions. The first one is the zeroth-order approximation of solvent effect,2-3 i.e.,

17

solvent is treated as continuum, and the solvent effect is reflected by the parameter setting of

18

dielectric constant during the calculation. In other words, the reaction is considered to occur in

19

a vacuum or gas phase or in a continuum. At present, most of the quantum chemistry

20

calculations are based on this approach.6-7 This method has the advantage of low computational

21

cost and can make reasonable explanation for gaseous phase reactions. However, in most

22

liquid phase reactions, especially when the solvent is composed of polar or charged molecules

23

and in addition the reactant and product molecules have polarity changes during the reaction 2


Page 2 of 35

Page 3 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

process, there is a significant and even qualitative difference in the calculated reaction free

2

energy with or without properly considering the solvent’s role.1, 8 This discrepancy is certainly

3

crucial to determine the reaction path and reaction mechanism.

4

The second category refers to the art of simulations coupling classical mechanics and

5

quantum mechanics (QM/MM), firstly proposed by Car and Parrinello, which takes the solvent

6

effect strictly into account through multiscale calculations.9 In principle, the reaction kinetics

7

including the reaction path can be exactly described by using QM/MM. The essential idea of

8

QM/MM method is to divide the reaction system into two regions.9-10 The core zone (quantum

9

system) undergoes an intrinsic reaction, which generally contains only the molecules and ions

10

involved in the reaction (sometimes including catalysts). The peripheral zone is a solvent

11

system that did not participate in the reaction, but interacts with the quantum system to affect

12

the occurrence of the reaction and the choice of the reaction path. Reaction substrates that have

13

not reacted yet are also included in the solvent system. The solvent system can be described by

14

molecular mechanics, whose influence on the quantum system is reflected in the crossover

15

Hamiltonian operator.10

16

In general, QM/MM calculation is theoretically complicated and computationally

17

expensive,11 and thus its extensive application has been limited. The zeroth-order

18

approximation may bring qualitative errors. Therefore, researchers introduce several more

19

advanced statistical mechanics theories to describe the solvent effect by combined quantum

20

mechanics and statistical mechanics.12-16 So far, a series of approaches have been developed,

21

including multi-configuration self-consistent field theory (MCSCF),17-18 the averaged solvent

22

electrostatic potential theory (ASEP),19-21 mean field approximation (MFA),22 and the 3



1

reference interacting site model self-consistent field (RISM/SCF) method.23 The most notable

2

method probably is the RISM/SCF theory which was developed based on RISM. RISM

3

represents an intermediate approach between the continuum description and molecular

4

simulation method, in which no solvent molecules are explicitly considered while the solvent

5

structure is captured by means of the solvent density distribution function.24-27 The solvent

6

density distribution can be obtained by solving the RISM equation iteratively coupled with a

7

closure relation. Within the framework of RISM/SCF, the solvent structure is affected by the

8

chemical nature and molecular configuration of solute (reaction subsystem) through the

9

solute-solvent interaction. Once the solute (namely, the chemical nature and/or molecular

10

configuration) changes during a reaction, the solvent structure is updated by solving the RISM

11

equation again. The system reaches the final self-consistent equilibrium between the reaction

12

sub-system and the solvent after multiple iterations.28

13

Whereas the RISM/SCF method can describe the solvent effect in a self-consistent

14

manner with reasonable computational cost,29 its extension to investigating the reaction

15

dynamics is restricted. This is because the involved RISM theory belongs to the scope of

16

integral equation methods, which are basically established through the correlation functions

17

of equilibrium systems. In other words, RISM is not tractable for describing dynamic

18

properties such as diffusion.

19

Similar to RISM theory, classical density functional theory (CDFT) represents an

20

efficient and robust modern statistical mechanics theory.30 Unlike quantum density functional

21

theory (QDFT), which relies on electron density distribution as functional variable, CDFT

22

takes the spatial density distribution of molecular species as functional variable, and it allows 4


Page 4 of 35

Page 5 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

us to predict the structure and thermodynamic properties of large-scale fluid systems involving

2

no chemical reaction.31 Depending on the geometrical complexity of system components,

3

three different versions of CDFT have been developed: atomic DFT, molecular DFT and

4

polymeric DFT.30 As reflected by their names, atomic DFT can be employed to investigate

5

simple fluid systems composed of spherical particles,22 while molecular DFT can be adopted

6

to investigate molecular fluid systems by taking into account of the orientational contribution

7

of non-spherical rigid molecules.32-33 In a further step, the chain-connectivity contribution in

8

polymeric molecules is considered in polymeric DFT.34-35 Although the functionals in these

9

three versions of CDFT are notably different, they share the same variational framework.30

10

Besides these sub-branches of CDFT, a special variant of molecular DFT, i.e., site DFT,

11

should be mentioned as well. Site DFT, originally proposed by Chandler et al.,36 can be

12

utilized to study molecular fluid systems by dealing with the molecules with their atomic site

13

densities.33, 37 However, limited by the constraint of rigid-molecular geometry, the site density

14

in site DFT is not an independent variable.33, 38

15

In recent decades, CDFT has been applied to tackle various problems ranging from

16

fundamental physics to engineering challenges. To name but a few, CDFT was extended for

17

accurate prediction of solvation free energy, and recently it has been adopted for fast

18

screening of hydration free energy.39-40 In parallel, the screening of gas adsorption in porous

19

materials was also conducted for targeting porous materials with high gas storage

20

capacity.41-42 The extension of CDFT to non-equilibrium field is also successful. The

21

dynamical version of atomic DFT, conventionally called dynamic DFT, has been extensively

22

developed and widely applied in the computational-aided design of supercapacitor with high

23

performance.43-46

5



1

In light of the success and robustness of CDFT, an interesting question arises on the

2

combination of CDFT with QDFT for assessing the solvent effects of chemical reactions. The

3

first attempt was made probably by Arias and his co-workers,47 who proposed a joint DFT

4

very recently. Within the framework of joint DFT for a solute-solvent system, the solvent is

5

commonly treated with the polarizable continuum model, in which the solvent density is

6

directly related with the solute electron density, rather than an independent variable.38, 48 In

7

addition, the communication between solute and solvent is described with the so-called

8

density-only electronic den7sity functional method,48-49 and this makes the electron density to

9

be the unique variable in the joint DFT, allowing for robust and efficient minimization

10

procedure.50 Although explicit description of water system was individually developed later

11

on along the scope of site DFT by the same group, the deviations of the predicted oxygen and

12

hydrogen density profiles from the experimental results are apparent,38 and furthermore the

13

incorporation of this extended site DFT into the joint DFT has been rarely reported.

14

Herein we report the first attempt on the combination of QDFT and molecular DFT in

15

different length scales to investigate the chemical reaction in aqueous solution. By using this

16

multiscale description, we propose a reaction density functional theory (RxDFT). To

17

demonstrate, the RxDFT is thereafter applied to describe the geometry and thermodynamics of

18

the tautomerization reaction of glycine in aqueous solution. The remainder of this work is

19

organized as follows: In next section, the theoretical framework of RxDFT is introduced,

20

along with the illustrative computational details for the tautomerization reaction. The results

21

are discussed in the sec. III, followed by a short conclusion in sec. IV.

22

II.

23

1.

Modelling and Theory Theoretical framework of RxDFT

6


Page 6 of 35

Page 7 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1 2

Figure 1. Thermodynamic circle for a representative chemical reaction (A+BC) occurs in

3

gas phase and in solution. FsolA , FsolB and FsolC are the solvation free energies of reagent and

4

product molecules, respectively, and E R denotes the reaction free energy in gas phase, and

5

FR is the reaction free energy in solution.

6 7

Considering a representative simple chemical reaction in a solvent: A+BC, where A

8

and B denote two different reagent molecules, and C denotes the product molecule. To

9

address the solvent effect, a thermodynamic circle can be established between the reaction in

10

dilute gas phase and that in the solvent, as depicted in Fig. 1. Since the solvation free energy

11

gives the reversible minimum work of moving a solute molecule from dilute gas phase to the

12

solvent,39 the solvation process of molecules A, B or C in the chemical reaction can be

13

C characterized with the individual solvation free energies, namely, FsolA , FsolB and Fsol . The

14

reaction free energy is denoted as FR in the solvent and E R in the gas phase. Note that in

15

dilute gas phase, the reaction free energy is usually approximated with the reaction energy as

16

the entropic contribution is overall neglectable. 7



1 2 3 4 5

Page 8 of 35

Following the thermodynamic circle, the reaction free energy in the solvent can be calculated by:

FR = ER + FsolC − FsolA − FsolB .

(1)

By setting ∆FS = FsolC − FsolA − FsolB , the above equation can be further simplified as:

FR = ER + ∆FS .

(2)

6

Here ∆FS represents the change in the solvation free energy after the occurrence of

7

chemical reaction. Clearly, when ∆FS is relatively small, eq.(2) recovers to the conventional

8

zeroth-order approximation of the solvent effect. However, if ∆FS is not ignorable, the

9

contribution of reaction free energy originates from two parts, as similar to the RISM/SCF

10

method,1 namely, the intrinsic reaction free energy and the solvation free energy.

11

Correspondingly, we may divide the system into two regions: the core region (quantum

12

system) undergoing an intrinsic reaction, and the peripheral region (solvent system). Whereas

13

several interesting attempts have been reported regarding the choice of the interface between

14

the core region and peripheral region,51-54 e.g., the solvent molecules in the first layer

15

surrounding the intrinsic reaction can be included into the core region; here we take the core

16

region as the one containing only the molecules and ions involved in the reaction.

17

The solvation process influences the mechanisms of chemical reactions in solution, but

18

the relaxation of the solvent is often much faster than the reactive event.4, 55-56 Namely,

19

solvation free energy calculation can be decoupled with the computation of the intrinsic

20

reaction free energy. While the intrinsic reaction free energy can be individually computed

21

with QDFT, the solvation free energy can be addressed with CDFT, being expressed as a

22

functional of the solvent density distribution function surrounding the solute: 8


Page 9 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

Fsol [ ρ ( x);Vext (x); µ , V , T ] = Θ[ ρ ( x);Vext (x); µ , V , T ] − Θ[ ρ 0 ;Vext (x); µ , V , T ] ,

2

where ρ (x) is the density profile of solvent molecules with x=(r, Ω) refers to the coordinates

3

and spatial orientation of solvent molecules. ρ 0 and µ are the solvent bulk density and bulk

4

chemical potential, respectively. V and T are the system volume and absolute temperature,

5

and they are given in prior and thus omitted below in the set of functional variables. Vext ( x )

6

represents the external potential to the solvent system, and it originates from the

7

solvent-solute interaction. Clearly, Vext ( x) changes when the molecular geometry and

8

charge distribution of the quantum system change duo to the occurrence of chemical reaction.

9

Certainly, Vext ( x ) is responsible for the interaction between the quantum system and solvent

10

system, and it comprises typically the coulomb interaction and Lennard-Jones interaction, as

11

disposed in other multiscale approaches.23, 47, 57 At a given Vext ( x ) , the grand potential can be

12

calculated as:

(3)

13

Θ[ ρ ( x);Vext ( x)] = k BT ∫ dx ρ ( x ){ln[ ρ ( x) Λ 3 ] − 1} + F exc [ ρ ( x )] + ∫ dx[Vext ( x) − µ ]ρ ( x ) ,

14

where k B is the Boltzmann constant; Λ is an effective thermal wavelength of solvent

15

molecule. F ex [ ρ ( x )] is the excess free energy functional accounting for the solvent

16

intermolecular interaction, of which no exact expression is available. To date, several

17

treatments have been developed for formulating the excess free energy functional. The

18

detailed formulation of this term will be discussed below.

(4)

19

The solvent density distribution minimizes the free energy functional of the solvent

20

system under the thermodynamic equilibrium condition. In other words, the solvent density

9



Page 10 of 35

1

distribution and the solvation free energy of the solute can be obtained by minimizing

2

solvation free energy functional in eq.(3),39, 58 namely:

3

δ Fsol [ ρ (x);Vext ( x)] =0 ⇒ δρ (x)

ρ ( x ) = ρ eq ( x) .

(5)

F = Fsol

4

Because CDFT directly relates the free energy to the density distribution of system

5

components and bypasses the redundant thermodynamic integral, it greatly reduces the

6

computational cost and displays great superiority compared with molecular simulations.58

7

With the help of CDFT, the hydration free energies of different neutral solutes and ions can

8

be quantitatively predicted in an efficient way.39 On the other hand, compared with the RISM

9

theory, the theoretical framework of CDFT is more flexible, and it has recently been extended

10 11

to non-equilibrium system to describe the diffusion dynamics. 2.

Reaction free energy and activation free energy

12

The evaluation of activation free energy and reaction free energy represents one of the

13

most important routes to explore reaction mechanisms. From the perspective of physical

14

chemistry, the reaction free energy supplies an important criterion for judging the pathway of

15

reaction, and the magnitude of the activation energy determines the reaction probability.

16

According to the transition state theory,59 as shown in Fig. 2a, the free energies of the quantum

17

system in gas phase before and after reaction are denoted as F0 and F2 , respectively.

18

Conventionally, the reaction free energy is calculated as ER = F2 − F0 . Assuming that the free

19

energy of transition state system is F1 , the energy barrier or activation free energy of reaction

20

is E A = F1 − F0 .

10


Page 11 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1 2

Figure 2. Typical energy profiles for a reaction a) in gas phase, and b) in a solvent. The

3

difference between the energies of transition state and reactant is the energy barrier or

4

activation free energy of reaction; it represents the minimum energy that a reacting system

5

must acquire for the transformation to take place. The difference between the energies of the

6

product and the reactant is the free energy of reaction.

7 8

Similarly, the reaction free energy is calculated as FR = F2′ − F0 in a solution, as shown

9

in the Fig. 2b. The energy barrier or activation free energy of reaction is FA = F1′ − F0 . Owning

10

to the decoupling between the intrinsic reaction free energy and solvent effect, FR

11

FA )

12

reaction with the solvation free energy difference, as stated in Eq.(2).

13

3. Revisit of the tautomerization reaction of glycine

(and

can be evaluated by summing up the corresponding E R (and E A ) for the intrinsic

14

Tautomerization is pervasive in chemistry and biology.60 There have been growing

15

theoretical interests in studying tautomerization in solution,57, 61-69 and this is because not

16

only

tautomerization reactions in solution can provide tractable model systems,70-72 but also

11



Page 12 of 35

1

these reactions are of practical importance, and they could be the rate-determining step for

2

kentonization and aldol condensation in biomass conversion and even decomposition of

3

molecular explosives.73-74 Herein we revisit of the tautomerization of glycine in water.

4

Glycine is the smallest amino acid that has been found in its neutral form (N) in the gas phase

5

or vacuum, while in aqueous solution the zwitterionic form (Z) dominates due to the

6

electrostatic interaction with the water molecules. Although the free energy change from the

7

neutral form to the zwitterionic form is about -7.3 kcal/mol measured in experiments,61, 75 there

8

are two possible reaction paths: the first path is the so called nonassisted concerted mechanism

9

in which the hydroxyl hydrogen (H6) forms a 5-membered ring transition state 1 (TS1) before

10

moving to the amino nitrogen (N1) of the zwitterionic form (Z) (see Fig. 3a). The second path

11

is the assisted concerted mechanism in which the glycine forms a 7-membered ring transition

12

state 2 (TS2) with water and then dehydrates to form a zwitterion structure (Z) (see Fig. 3b).

a)

H6 N1

O4 C3

C2

O5

N

TS1

Z

TS2

Z

Ow

b) Hw H6 N1

O4 C3

C2

O5

13

N

14

Figure 3. Schematic representations and atomic numbering for the neutral, transition, and

15

zwitterionic structure of glycine with: a) nonassisted concerted mechanism and b) assisted

16

concerted mechanism. 12


Page 13 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1 2

The free energy profile for the Glycine tautomerization reaction in aqueous solution is

3

depicted in Fig. 4. In the initial state, the grand potential of the system is composed the intrinsic

4

free energy of the quantum system and grand potential energy of the surrounding solvent

5

system, with the external potential corresponding to the molecular geometry and charge

6

distribution in the quantum system. Similar composition holds regarding the grand potentials

7

for the transition state system and for the final state. In each state, the molecular structure and

8

charge distribution of glycine are determined upon quantum DFT calculation, which hereby

9

provide the necessary input information, through the solute-solvent interaction potential, for

10

the CDFT calculation of the hydration free energies.

11 12

Figure 4. Schematic free energy profile for glycine tautomerization reaction in aqueous

13

solution.

14 15

III.

16

1. Intrinsic free energy calculation with QDFT

Calculation Methods and Details

13



Page 14 of 35

1

To calculate intrinsic reaction free energy ( E R ), we first optimize the structure of

2

glycine in water by treating the solution with polarizable continuum model, and then the

3

obtained geometric configuration is utilized for a single-point frequency calculation in gas

4

phase for determining E R . The quantum calculations are carried out by using the

5

GAUSSIAN 09 package of programs76 at the MP2,77 B3LYP78 and M06-2X79 levels,

6

respectively. Following previous studies of glycine,57, 64, 69 the 6-31+G** basis set has been

7

chosen, and this basis set includes the diffuse functions on all atoms except the hydrogens

8

and the polarization functions on all atoms. Molecular geometries of the reactants, products

9

and transition states (TS1 and TS2) have been optimized using Berny’s algorithms.80 A

10

single-point frequency calculation is followed to ensure that the final structure is obtained

11

without imaginary frequency (reactant and product) or with only a single imaginary

12

frequency (transition states). These TS1 and TS2 transition states are verified by performing a

13

normal-mode analysis, which provides a single imaginary frequency ( υTS1 = -1179 cm-1; υTS2

14

= -1313 cm-1 for the glycine at MP2 level). Intrinsic reaction coordinate (IRC)81 calculations

15

are performed for ensuring the transition states connect to the correct reactants and products.

16

2. Hydration free energy calculation with CDFT

17

2.1 excess free-energy functional

18 19

Following our previous work39, the formulation of the excess free energy functional for inhomogeneous water is formulated as:

20

1 ∆F exc [ ρ (x)]= − k BT ∫∫ ∆ρ ( x1 ) c(2) (x1 , x 2 ;[ ρ0 ])∆ρ ( x 2 ) dx1dx 2 + F B [ ρ ( x )] , 2

21

where ∆ρ (x) = ρ (x) − ρ 0 . On the right hand side of eq.(6), the first term represents the

22

homogeneous reference fluid (HRF) approximation. It amounts to a second-order Taylor

23

expansion of the excess free-energy functional for inhomogeneous water around the bulk one,

14


(6)

Page 15 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

which is known to be equivalent to the hypernetted-chain (HNC) approximation in integral

2

(2) equation theories.22 c (x1 , x2 ;[ρ0 ]) is the two-body direct correlation function (DCF) of the

3

bulk water under the same thermodynamic condition, and it can be extracted from molecular

4

dynamics simulation.82 In this work, the DCF of bulk water serves as an input of the density

5

functional, prepared elaborately in prior.83 The second term in eq.(6) represents the

6

contribution from the multi-body correlations, conventionally designated as bridge functional,

7

which can be rigorously expressed as a series of expansions in terms of the three-body and

8

higher-order correlations of the bulk water. As in many work,39, 82, 84 we here approximate the

9

bridge functional with the one for a reference inhomogeneous hard-sphere fluid system with

10

density distribution ρ (r) :

F B [ ρ ( x )] ≈ FhsB [ ρ ( r )]

11

(7)

12

Here ρ (r ) = ∫ ρ ( x)d Ω . Similar to earlier work,85 the bridge functional of the corresponding

13

hard sphere system is calculated from the modified fundamental measure theory (MFMT)86: exc exc FhsB [ ρ ( r )]=Fhsexc [ ρ ( r )] − Fhs-b ( ρ b ) − µ hs-b ∫ ∆ρ ( r ) d r

14

1 (2) + k BT ∫∫ ∆ρ ( r1 ) chs-b (| r1 − r2 |) ∆ρ ( r2 ) dr1dr2 2

(8)

15

where Fhsexc [ ρ ( r )] is the excess free-energy functional of the inhomogeneous hard-sphere

16

exc exc (2) system; Fhs-b , µ hs-b and chs-b are the excess free-energy, excess chemical potential and

17

(2) DCF for the corresponding hard-sphere bulk system. In this work, Fhsexc [ ρ ( r )] and chs-b are

18

exc exc calculated from MFMT,86 while Fhs-b and µ hs-b are calculated from the Carnahan-Starling

19

equation of state.87 The effective hard-sphere diameter of the reference system is 2.95Å for

20

SPC/E water.33

21

2.2 Numerical aspects

22

The minimization of eq.(5) gives rise to the Euler-Lagrange equation: 15



1

 ρ (x)  δ∆F exc [ ρ ( x)] k BT ln  = 0.  + Vext ( x) + δρ ( x)  ρ0 

Page 16 of 35

(9)

2

The inhomogeneous density ρ (x) is projected onto an orthorhombic position grid of

3

N x × N y × N z nodes with periodic boundary conditions applied in all three directions. The

4

variational density ρ (x) , subject to the minimization of F [ ρ (x);Vext (x)] , is optimized by the

5

limited memory Broyden-Fletcher-Goldfarb-Shanno (L-BFGS) method.88 This quasi-Newton

6

algorithm requires the formulation of the functional F [ ρ (x);Vext (x)] and its first derivative

7

with respect to the density at each node. The convolutions in F exc [ ρ (x)] are calculated in

8

k-space by using fast Fourier transforms (FFT) as implemented with the FFTW3 library.89

9

The temperature T = 298.15 K, and water bulk density 0.0333 particles/Å3, corresponding to

10

the mass concentration 0.997 g/cm3.

11

The external potential Vext (x) , involving atomic Lennard-Jones (LJ) and partial charges

12

parameters, is computed only once before the numerical minimization. In each calculation, a

13

single glycine molecule (i.e. in neutral, zwitterionic or transition state form) is fixed at the

14

center of a cubic grid box with box length of approximately 40 Å, and the external potential

15

to the water molecule locating at r with orientation Ω sums up the pair interactions

16

between the site i in glycine and site j in water, reading:

17

 σ Vext ( r, Ω ) = ∑∑ 4ε ij  ij  rij i j  a

b

12 6   σ ij   qi q j .  −    +   rij   4πε 0 rij

(10)

18

Here rij is the interatomic distance between site i and site j , while qi and q j are the

19

partial charges on these sites. The conventional Lorentz-Berthelot combination rule is used to

20

generate the LJ parameters ε ij and σ ij for the unlike pairs.

16


Page 17 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

For water molecule, we have used the SPC/E model potential parameters.90 The charge

2

distributions in glycine are calculated by fitting the individual (atomic) charges to the natural

3

population analysis (NPA) or molecular electrostatic potential. The latter are computed from

4

the QM wave function at series of points selected by the charges from electrostatic potential

5

(CHELPG) algorithm.91 The LJ parameters of the nonassisted and assisted mechanisms for

6

neutral, zwitterionic and transition state glycine are taken from the OPLS force field database

7

(Tables S1-S6).92

8

IV.

9

1. Molecular Geometric Configurations

Results and Discussion

10

For the nonassisted concerted mechanism, the neutral (N), zwitterionic (Z) and transition

11

state 1 (TS1) conformations of glycine in water are displayed in Fig.3a. In general, the neutral

12

glycine presents a large number of possible conformers in various environments,62, 93-95 and its

13

conformer in water is determined by minimizing the absolute value of its intrinsic free energy

14

by treating water as polarizable continuum.15 Note that the obtained structure of glycine in

15

water is different to its counterpart in gas phase. The same method is applied to determine the

16

structures of glycine in zwitterionic and TS1 states in water. In gas phase, the geometric

17

configurations of zwitterionic and TS1 glycine molecules can’t be located through the

18

minimization of intrinsic free energy,63 and thus the neutral, zwitterionic and transition state 1

19

conformations of glycine are assumed to be identical with their counterparts in water.

20

Fig. 3b presents the configuration involving only one water molecule in the assisted

21

concerted mechanism. The hydroxyl group is oriented to the amino group, which facilitates the

22

formation of intramolecular hydrogen bonds by promoting the proton transfer of H6 atom to

23

Ow atom (assisted mechanism) or N1 atom (nonassisted mechanism). The geometric

24

configurations of glycine in neutral, zwitterionic and transition state 2 (TS2) in the assisted 17



1

concerted mechanism are determined with the similar method as in the nonassisted concerted

2

mechanism.

3

The geometric parameters described for the neutral (N), zwitterionic (Z) and transition

4

state (TS1 and TS2) of glycine are given in Table 1. All these optimized structures are

5

calculated in the same computational level mentioned in section III. These results are

6

consistent with the observations by other workers for amino acids in water.64-67

7

18


Page 18 of 35

Page 19 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47


Table 1. Geometric parameters for glycine at the different computational levels used.

N C2-N1 C2-C3 C3-O5 C3-O4 C2-H7 O4-H6 N1-H6 O4-N1 H6-Ow Ow-Hw Hw-N1 N1-C2-C3 O5-C3-C2 O4-C3-C2 N1-H6-O4 O4-H6-Ow H6-Ow-Hw Ow-Hw-N1 N1-C2-C3-O4 C3-C2-N1-H6 C3-C2-N1-Hw

TS1

MP2

B3LYP

M06-2X

MP2

1.46 1.53 1.23 1.34 1.09 1.00 2.57 109.8 123.6 113.8 128.8 0.0 0.0 -

1.47 1.53 1.22 1.33 1.09 1.00 2.57 109.8 123.9 113.4 128.8 0.0 0.0 -

1.46 1.53 1.21 1.33 1.09 0.99 2.57 109.9 123.5 113.9 127.3 0.0 0.0 -

1.48 1.54 1.24 1.31 1.09 1.24 1.27 2.37 103.8 121.8 110.6 140.8 0.0 0.0 -

TS2

B3LYP M06-2X

1.49 1.55 1.23 1.30 1.09 1.26 1.27 2.37 103.8 121.8 110.6 140.0 0.0 0.0 -

1.48 1.55 1.23 1.29 1.09 1.30 1.22 2.37 103.5 120.9 110.9 139.3 0.0 0.0 -

MP2

1.48 1.53 1.24 1.31 1.09 1.17 2.89 1.25 1.24 1.25 114.1 117.6 118.9 170.3 86.0 162.3 -36.0 52.2

19


Z

B3LYP M06-2X

1.49 1.54 1.23 1.30 1.09 1.19 2.90 1.22 1.24 1.26 115.2 117.1 119.4 169.8 86.7 162.3 -32.2 47.8

1.48 1.54 1.23 1.29 1.09 1.24 2.87 1.17 1.22 1.27 114.4 117 118.8 168.1 87.9 161.2 -32.6 51.3

MP2

1.50 1.55 1.26 1.27 1.09 1.04 2.56 107.7 115.5 115.1 123.3 -0.5 0.7 -

B3LYP M06-2X

1.50 1.55 1.25 1.27 1.09 1.05 2.55 107.7 116.0 114.9 124.2 -0.3 0.4 -

1.49 1.55 1.24 1.26 1.09 1.05 2.55 107.7 115.6 115.0 122.7 0.0 0.0 -


1

In relation to the TS1 structure on the way from neutral to zwitterionic, O4-H6 distance

2

is much larger at M06-2X level (0.31 Å) than at the MP2 (0.24 Å) or B3LYP (0.26 Å). In all

3

cases, the process puts the H6 at a distance of about 1.2 Å of the N1 by extending the O4-H6

4

bond. Moreover, we find that the TS1 glycine is a nearly planar 5-membered ring with very

5

similar distances and bond angles.

6

When a water molecule is involved as part of the reaction mechanism, we find that the

7

Hw links to N1, and the H6 links to Ow. We also observe in Table 1 that the length of the

8

Ow-Hw bond increases to 1.24 Å (MP2 and B3LYP) or 1.22 Å (M06-2X), with the Ow-H6

9

bond length at 1.25 Å (MP2), 1.22 Å (B3LYP) and 1.17 Å (M06-2X), and the O4-H6

10

distance is slightly increased, which allows to release the H6-Ow-Hw water molecule. In the

11

TS2 structure, the formation of a 7-membered nonplanar ring renders the structure less

12

stressed than the TS1 structure, presenting a larger bond angle especially the O4-H6-Ow and

13

Ow-Hw-N1 angles involved in the proton transfer.

14

2. Energetics and Reaction Path

15

Table 2 and Fig. 5 present the relative energies and free energy profiles for glycine at

16

the neutral, zwitterionic, and transition states in water, with different levels of QDFT

17

calculations. In the gas phase for both the assisted and nonassisted mechanisms, the

18

tautomerization process N→Z is clearly unfavorable, with higher activation barriers and very

19

endothermic reaction energies. In aqueous solution, there is a significant reduction of the

20

activation barrier, generating depressed reaction energies of the process being obtained (they

21

even are exothermic for both mechanisms). This is because the dipole moment of glycine

22

increases along the reaction coordinate,72 leading to a continuous decrease of the hydration free

23

energy, and this suggests that the solvent plays an important role in stabilizing the final state.

24

Furthermore, the reduction of the activation barriers in the water-assisted mechanism is

25

slightly more pronounced than that in the non-assisted mechanism, indicating that the

20 ACS Paragon Plus Environment

Page 20 of 35

Page 21 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

solvent’s role is more significant in the water-assisted mechanism. Indeed, the dipole moment

2

of TS2 is 8.4 D, being larger than 8.2 D for TS1 at the B3LYP level. Tolosa et al.57 proposed a

3

scheme to describe the chemical reaction in water by means of molecular dynamics (MD)

4

simulation combining with quantum calculation at MP2 level, and the MD/MP2 method was

5

applied to study the tautomerization of four amino acids in water using free energy curves at

6

298 K. In their calculation, the TIP3P water model was employed.96

7

The values of the activation and reaction free energies with the nonassisted mechanism

8

and assisted mechanism predicted by RxDFT are quantitatively comparable with the

9

calculation results by Tolosa et al.57 Moreover, the combination of B3LYP with MDFT in the

10

nonassisted mechanism gives very satisfactory predictions (7.83 kcal/mol for the activation

11

free energy and -4.84 kcal/mol for the reaction free energy), compared with Tolosa’s results

12

(6.64 kcal/mol and -4.91 kcal/mol). Besides, the theoretical predictions are overall consistent

13

with the experimental measurements.68, 75, 61, 97

14

Benchmarked with the results from Tolosa, the free energy profiles of glycine in gas

15

phase can be satisfactorily predicted with these three methods, rationalizing our treatments on

16

the geometric structures of glycine in gas phase. In the aqueous solution, the results from

17

B3LYP + MDFT and M06-2X + MDFT method are better than that from MP2 + MDFT

18

method for the nonassisted mechanism, while the free energy profile from MP2 + MDFT

19

method is better than those from B3LYP + MDFT and M06-2X + MDFT methods for the

20

assisted mechanism.

21 22

Table 2. Relative energies (in kcal/mol) for the structures of glycine studied at the different

23

computational levels used. N nonassisted

gas

MP257

TS

0.00 11.56


Z 20.65


Page 22 of 35

MP2

0.00

9.04

20.47

B3LYP

0.00

8.61

20.73

M06-2X

0.00 10.17

23.53

MD/MP257

0.00

6.64

-4.91

MP2 + MDFT

0.00

6.24

-11.38

B3LYP + MDFT

0.00

7.83

-4.84

M06-2X + MDFT

0.00

9.62

-3.08

MP257

0.00 13.40

13.23

MP2

0.00 11.26

11.93

B3LYP

0.00

9.25

14.00

M06-2X

0.00

9.65

15.39

MD/MP257

0.00

5.29

-4.27

MP2 + MDFT

0.00

4.04

-5.85

B3LYP + MDFT

0.00

3.93

-9.08

M06-2X + MDFT

0.00

3.18

-7.22

0.00

6.90

-7.30

solution

gas

assisted

solution

Experimental 1

2 3

Figure 5. Free energy profiles for glycine using a) a nonassisted mechanism, b) an assisted


Page 23 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

mechanism. The label G denotes the reaction in gas phase and S denotes the reaction in

2

aqueous solution. The free energy profiles of MD/MP2 method are taken from Ref. 57.

3 4

V.

Conclusions

5

We have proposed an RxDFT by combining QDFT and CDFT, which enable us to

6

address the solvent effect of chemical reactions explicitly. For demonstration, the RxDFT is

7

applied to investigate the reaction paths associated with the tautomerization reaction of

8

glycine in water. In particular, the intrinsic reaction free energy and solvent effect are

9

decoupled, and thereafter treated separately by using QDFT and CDFT. The geometric

10

configurations and intrinsic free energies of and the charge distributions in glycine at

11

different reaction states have been computed with QDFT at the MP2/6-31+G**,

12

B3LYP/6-31+G** and M06-2x/6-31+G** computational level. The obtained geometries and

13

charge distributions determine the glycine-water interaction, which serves, hereby, as an input

14

for the CDFT calculation for determining the change in solvation free energy during the

15

occurrence of tautomerization. The summation of the intrinsic reaction free energy with the

16

solvation free energy difference gives the total reaction free energy in water.

17

The intrinsic free energies are obtained at the MP2, B3LYP and M06-2x computational

18

levels, and the profile confirms that the chemical reaction in gas phase is endothermic. The

19

final free-energy curves for both nonassisted concerted mechanism and assisted concerted

20

mechanism are obtained by using RxDFT, showing that the reaction in water is exothermic.

21

The theoretical predictions from RxDFT are compared with the reported MD/MP2

22

calculations and experimental measurement, and an overall satisfactory agreement has been

23

found, demonstrating that the proposed RxDFT is accurate and reliable.

24

Although water is uniquely considered as a common representative for reaction

25

solvents in this work, the proposed protocol in RxDFT can be extended to address the



1

chemical reactions in other solvents. Due to the low computational cost, we consider RxDFT

2

may provide an efficient alternative to ab initio correlated methods and QM/MM methods for

3

accurate theoretical study of reaction free energies in solvents.

4 5

Associated Content

6

Supporting information

7

The atomic (i.e. site) coordinates, Lennard-Jones parameters and charges for the glycine

8

molecules used in the RxDFT calculation.

9 10

Acknowledgements

11

This work is supported by National Natural Science Foundation of China (No. 91434110 and

12

21878078), National Natural Science Foundation of China for Innovative Research Groups

13

(No. 51621002), and the 111 Project of China (No. B08021). SZ acknowledges the support of

14

Fok Ying Tong Education Foundation (151069).

15 16

References

17

1.

18 19

24. 2.

20 21

3.

Tomasi, J.; Mennucci, B.; Cammi, R., Quantum Mechanical Continuum Solvation Models. Chem. Rev. 2005, 105, 2999-3093.

4.

24 25

Mennucci, B.; Cammi, R., Continuum Solvation Models in Chemical Physics: From Theory to Applications; John Wiley & Sons, 2007.

22 23

Hirata, F., Molecular Theory of Solvation; Springer Science & Business Media, 2003; Vol.

Reichardt, C.; Welton, T., Solvents and Solvent Effects in Organic Chemistry; John Wiley & Sons, 2011.

5.

Carey, F. A., Organic Chemistry; The McGraw-Hill Companies, 2004.


Page 24 of 35

Page 25 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

6.

2 3

Miertus,̃ S.; Tomasi, J., Approximate Evaluations of the Electrostatic Free Energy and Internal Energy Changes in Solution Processes. Chem. Phys. 1982, 65, 239-245.

7.

Miertuš, S.; Scrocco, E.; Tomasi, J., Electrostatic Interaction of a Solute with a Continuum.

4

A Direct Utilizaion of Ab Initio Molecular Potentials for the Prevision of Solvent Effects.

5

Chem. Phys. 1981, 55, 117-129.

6

8.

Prototype SN2 Reaction in Solution: Cl− Attack on CH3Cl. J. Chem. Phys. 2014, 140, 289.

7 8 9

Kuechler, E. R.; York, D. M., Quantum Mechanical Study of Solvent Effects in a

9.

Car,

R.;

Parrinello,

M.,

Unified

Approach

for

Molecular

Dynamics

and

Density-Functional Theory. Phys. Rev. Lett. 1985, 55, 2471.

10

10. Vreven, T.; Frisch, M. J.; Kudin, K. N.; Schlegel, H. B.; Morokuma, K., Geometry

11

Optimization with QM/MM Methods II: Explicit Quadratic Coupling. Mol. Phys. 2006,

12

104, 701-714.

13 14 15 16

11. Teodoro, L.; Fawzi, M.; Alessandro, L.; Parrinello, M., An Efficient Real Space Multigrid QM/MM Electrostatic Coupling. J. Chem. Theory Comput. 2005, 1, 1176-1184. 12. Orozco, M.; Luque, F. J., Theoretical Methods for the Description of the Solvent Effect in Biomolecular Systems. Chem. Rev. 2000, 100, 4187-4226.

17

13. Muñoz, J.; Barril, X.; Luque, F. J.; Gelpí, J. L.; Orozco, M., Partitioning of Free Energies

18

of Solvation into Fragment Contributions: Applications in Drug Design; Springer US,

19

2001, p 143-168.

20 21

14. Vreven, T.; Morokuma, K., On the Application of the IMOMO (Integrated Molecular Orbital + Molecular Orbital) Method. J. Comput. Chem. 2015, 21, 1419-1432.

22

15. Vreven, T.; Byun, K. S.; Komáromi, I.; Dapprich, S.; Montgomery, J. A.; Morokuma, K.;

23

Frisch, M. J., Combining Quantum Mechanics Methods with Molecular Mechanics

24

Methods in ONIOM. J. Chem. Theory Comput. 2006, 2, 815-826.



1

16. Mikkelsen, K. V.; Jo/Rgensen, P.; Jensen, H. J. R. A., A Multiconfiguration

2

Self-Consistent Reaction Field Response Method. J. Chem. Phys. 1994, 100, 6597-6607.

3

17. Bishop, D. M., Aspects of Non-Linear-Optical Calculations. Adv. Quantum Chem. 1994,

4 5 6

25, 1-45. 18. Shelton, D. P.; Rice, J. E., Measurements and Calculations of the Hyperpolarizabilities of Atoms and Small Molecules in the Gas Phase. Chem. Rev. 1994, 94, 3-29.

7

19. Galván, I. F.; Sánchez, M. L.; Martı́N, M. E.; Valle, F. J. O. D.; Aguilar, M. A., ASEP/MD:

8

A Program for the Calculation of Solvent Effects Combining QM/MM Methods and the

9

Mean Field Approximation. Comput. Phys. Commun. 2003, 155, 244-259.

10

20. Sanchez, M. L.; Aguilar, M. A.; Olivares, d. V., F. J, Study of Solvent Effects by Means of

11

Averaged Solvent Electrostatic Potentials Obtained from Molecular Dynamics Data. J.

12

Comput. Chem. 2015, 18, 313-322.

13

21. Sánchez, M. L.; Martín, M. E.; Aguilar, M. A.; Olivares, d. V., F. J, Solvent Effects by

14

Means of Averaged Solvent Electrostatic Potentials: Coupled Method. J. Comput. Chem.

15

2015, 21, 705-715.

16

22. Hansen, J. P.; Mcdonald, I. R., Theory of Simple Liquids (Second Edition), 1990.

17

23. Hirata, F.; Sato, H.; Ten-No, S.; Kato, S., The RISM-SCF/MCSCF Approach for

18

Chemical Processes in Solutions. Eastern Hemisphere Distribution 2001, 417.

19

24. Andersen, H. C.; Chandler, D., Optimized Cluster Expansions for Classical Fluids. I.

20

General Theory and Variational Formulation of the Mean Spherical Model and Hard

21

Sphere Percus-Yevick Equations. J. Chem. Phys. 1972, 57, 1918-1929.

22 23

25. Chandler, D.; Andersen, H. C., Optimized Cluster Expansions for Classical Fluids. II. Theory of Molecular Liquids. J. Chem. Phys. 1972, 57, 1930-1937.


Page 26 of 35

Page 27 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

26. Andersen, H. C.; Chandler, D.; Weeks, J. D., Optimized Cluster Expansions for Classical

2

Fluids. III. Applications to Ionic Solutions and Simple Liquids. J. Chem. Phys. 1972, 57,

3

2626-2631.

4 5

27. Hudson, S.; Andersen, H. C., Optimized Cluster Expansions for Classical Fluids. IV. Primitive Model Electrolyte Solutions. J. Chem. Phys. 1974, 60, 2188-2188.

6

28. Wesolowski, T. A.; Weber, J., Kohn-Sham Equations with Constrained Electron Density:

7

An Iterative Evaluation of the Ground-State Electron Density of Interacting Molecules.

8

Chem. Phys. Lett. 1996, 248, 71-76.

9

29. Naka, K.; Sato, H.; Morita, A.; Hirata, F.; Kato, S., Rism-Scf Study of the Free-Energy

10

Profile of the Menshutkin-Type Reaction NH3 + CH3Cl → NH3CH3+ + Cl− in Aqueous

11

Solution. Theor. Chem. Acc. 1999, 102, 165-169.

12

30. Zhao, S.; Liu, Y.; Chen, X.; Lu, Y.; Liu, H.; Hu, Y., Chapter One-Unified Framework of

13

Multiscale Density Functional Theories and Its Recent Applications. Advances in

14

Chemical Engineering 2015, 47, 1-83.

15 16 17 18

31. Wu, J. Z., Density Functional Theory for Liquid Structure and Thermodynamics; Springer Berlin Heidelberg, 2009. 32. Zhao, S. L.; Liu, Y.; Liu, H. L.; Wu, J. Z., Site-Site Direct Correlation Functions for Three Popular Molecular Models of Liquid Water. J. Chem. Phys. 2013, 139, 064509.

19

33. Liu, Y.; Zhao, S.; Wu, J., A Site Density Functional Theory for Water: Application to

20

Solvation of Amino Acid Side Chains. J. Chem. Theory Comput. 2013, 9, 1896-1908.

21

34. Lian, C.; Cai, C.; Shen, X.; Zhao, S.; Yu, X.; Liu, H., Improved Oxidation of Hydrogen

22

Off-Gas by Hydrophobic Surface Modification: A Multiscale Density Functional Theory

23

Study. Particuology 2018.

24 25

35. Wu, J., Classical Density Functional Theory for Molecular Systems; Springer Singapore, 2017.



1

36. Chandler, D.; McCoy, J. D.; Singer, S. J., Density Functional Theory of Nonuniform

2

Polyatomic Systems. I. General Formulation. J. Chem. Phys. 1986, 85, 5971-5976.

3

37. Reddy, G.; Lawrence, C. P.; Skinner, J. L.; Yethiraj, A., Liquid State Theories for the

4

Structure of Water. J. Chem. Phys. 2003, 119, 13012-13016.

5

38. Sundararaman, R.; Arias, T. A., Efficient Classical Density-Functional Theories of

6

Rigid-Molecular Fluids and a Simplified Free Energy Functional for Liquid Water.

7

Comput. Phys. Commun. 2014, 185, 818-825.

8 9 10 11

39. Zhao, S.; Jin, Z.; Wu, J., New Theoretical Method for Rapid Prediction of Solvation Free Energy in Water. J. Phys. Chem. B 2011, 115, 6971-6975. 40. Wu, H.; Li, Y.; Kadirov, D.; Zhao, S.; Lu, X.; Liu, H., An Efficient Molecular Approach for Quantifying Solvent-Mediated Interaction. Langmuir 2017, 33, 11817-11824.

12

41. Liu, Y.; Liu, H.; Hu, Y.; Jiang, J., Development of a Density Functional Theory in

13

Three-Dimensional Nanoconfined Space: H2 Storage in Metal-Organic Frameworks. J.

14

Phys. Chem. B 2009, 113, 12326-12331.

15 16

42. Liu, Y.; Liu, H.; Hu, Y.; Jiang, J., Density Functional Theory for Adsorption of Gas Mixtures in Metal-Organic Frameworks. J. Phys. Chem. B 2010, 114, 2820.

17

43. Cheng, L.; Zhao, S.; Liu, H.; Wu, J., Time-Dependent Density Functional Theory for the

18

Charging Kinetics of Electric Double Layer Containing Room-Temperature Ionic Liquids.

19

J. Chem. Phys. 2016, 145, 204707.

20 21 22 23

44. Zhao, S.; Wu, J., Self-Consistent Equations Governing the Dynamics of Nonequilibrium Colloidal Systems. J. Chem. Phys. 2011, 134, 054514. 45. Ma, M.; Zhao, S.; Liu, H.; Xu, Z., Microscopic Insights into the Efficiency of Capacitive Mixing Process. AlChE J. 2017, 63, 1785-1791.


Page 28 of 35

Page 29 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

46. Ma, M.; Zhao, S.; Xu, Z., Investigation of Dielectric Decrement and Correlation Effects

2

on Electric Double-Layer Capacitance by Self-Consistent Field Model. Communications

3

in Computational Physics 2016, 20, 441-458.

4 5 6 7 8 9

47. Petrosyan, S. A.; Briere, J. F.; Roundy, D.; Arias, T. A., Joint Density-Functional Theory for Electronic Structure of Solvated Systems. Physical Review B 2007, 75. 48. Sundararaman, R. Joint Density-Functional Methods for First-Principles Chemistry in Solution. Cornell University, 2013. 49. Letchworth-Weaver, K.; Sundararaman, R.; A. Arias, T., First Principles Free-Energy Theory of Solvation with Atomic Scale Liquid Structure, 2017.

10

50. Sundararaman, R.; Letchworth-Weaver, K.; Schwarz, K. A.; Gunceler, D.; Ozhabes, Y.;

11

Arias, T. A., JDFTx: Software for Joint Density-Functional Theory. SoftwareX 2017, 6,

12

278-284.

13

51. Deng, S.; Cai, W., Extending the Fast Multipole Method for Charges inside a Dielectric

14

Sphere in an Ionic Solvent: High Order Image Approximations for Reaction Fields.

15

Journal of Computational Physics 2007, 227, 1246-1266.

16

52. Lin, Y.; Baumketner, A.; Deng, S.; Xu, Z.; Jacobs, D.; Cai, W., An Image-Based Reaction

17

Field Method for Electrostatic Interactions in Molecular Dynamics Simulations of

18

Aqueous Solutions. J. Chem. Phys. 2009, 131, 10B608.

19 20

53. Lin, Y.; Baumketner, A.; Song, W.; Deng, S.; Jacobs, D.; Cai, W., Ionic Solvation Studied by Image-Charge Reaction Field Method. J. Chem. Phys. 2011, 134, 01B627.

21

54. Qin, P.; Xu, Z.; Cai, W.; Jacobs, D., Image Charge Methods for a Three-Dielectric-Layer

22

Hybrid Solvation Model of Biomolecules. Communications in Computational Physics

23

2009, 6, 955-977.



1

55. Berg, M.; Harris, A. L.; Harris, C. B., Rapid Solvent-Induced Recombination and Slow

2

Energy Relaxation in a Simple Chemical Reaction: Picosecond Studies of Iodine

3

Photodissociation in CCl4. Phys. Rev. Lett. 1985, 54, 951-954.

4

56. Rosenthal, S. J.; Xie, X.; Du, M.; Fleming, G. R., Femtosecond Solvation Dynamics in

5

Acetonitrile: Observation of the Inertial Contribution to the Solvent Response. J. Chem.

6

Phys. 1991, 95, 4715-4718.

7

57. Tolosa, S.; Hidalgo, A.; Sansón, J. A., Amino Acid Tautomerization Reactions in

8

Aqueous Solution Via Concerted and Assisted Mechanisms Using Free Energy Curves

9

from Md Simulation. J. Phys. Chem. B 2012, 116, 13033-13044.

10 11 12 13 14 15

58. Zhao, S.; Ramirez, R.; Vuilleumier, R.; Borgis, D., Molecular Density Functional Theory of Solvation: From Polar Solvents to Water. J. Chem. Phys. 2011, 134, 194102. 59. Laidler, K. J.; King, M. C., Development of Transition-State Theory. J. Phys. Chem. 1983, 87, 2657-2664. 60. Smith, M. B.; March, J., March's Advanced Organic Chemistry: Reactions, Mechanisms, and Structure; John Wiley & Sons, 2007.

16

61. Wada, G.; Tamura, E.; Okina, M.; Nakamura, M., On the Ratio of Zwitterion Form to

17

Uncharged Form of Glycine at Equilibrium in Various Aqueous Media. Bull. Chem. Soc.

18

Jpn. 1982, 55, 3064-3067.

19 20

62. Császár, A. G., On the Structures of Free Glycine and α-Alanine. J. Mol. Struct. 1995, 346, 141-152.

21

63. Tortonda, F. R.; Pascual-Ahuir, J. L.; Silla, E.; Tuñón, I., Why Is Glycine a Zwitterion in

22

Aqueous Solution? A Theoretical Study of Solvent Stabilising Factors. Chem. Phys. Lett.

23

1996, 260, 21-26.


Page 30 of 35

Page 31 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

64. Tortonda, F. R.; Pascualahuir, J. L.; Silla, E.; Tuñón, I.; Ramı́Rez, F. J., Aminoacid

2

Zwitterions in Solution: Geometric, Energetic, and Vibrational Analysis Using Density

3

Functional Theory-Continuum Model Calculations. J. Chem. Phys. 1998, 109, 592-603.

4

65. Balta, B.; Aviyente, V., Solvent Effects on Glycine. I. A Supermolecule Modeling of

5

Tautomerization Via Intramolecular Proton Transfer. J. Comput. Chem. 2003, 24,

6

1789-802.

7

66. Tortonda, F. R.; Silla, E.; Tuñón, I.; Rinaldi, D.; Ruiz-López, M. F., Intramolecular Proton

8

Transfer of Serine in Aqueous Solution. Mechanism and Energetics. Theor. Chem. Acc.

9

2000, 104, 89-95.

10

67. Bowman, J. C.; Shadwick, B. A.; Morrison, P. J., Theoretical Study of the

11

Tautomeric/Conformational Equilibrium of Aspartic Acid Zwitterions in Aqueous

12

Solution. J. Phys. Chem. A 2000, 104, 6834-6843.

13 14 15 16

68. Slifkin, M.; Ali, S., Thermodynamic Parameters of the Activation of Glycine Zwitterion Protonation Reactions. J. Mol. Liq. 1984, 28, 215-221. 69. Tuñón, I.; Silla, E.; Ruiz-López, M. F., On the Tautomerization Process of Glycine in Aqueous Solution. Chem. Phys. Lett. 2000, 321, 433-437.

17

70. Senn, H. M.; Margl, P. M.; Schmid, R.; Ziegler, T.; Blöchl, P. E., Ab Initio Molecular

18

Dynamics with a Continuum Solvation Model. J. Chem. Phys. 2003, 118, 1089-1100.

19

71. Bandyopadhyay, P.; Gordon, M. S.; Mennucci, B.; Tomasi, J., An Integrated Effective

20

Fragment-Polarizable Continuum Approach to Solvation: Theory and Application to

21

Glycine. J. Chem. Phys. 2002, 116, 5023-5032.

22 23

72. Leung, K.; Rempe, S. B., Ab Initio Molecular Dynamics Study of Glycine Intramolecular Proton Transfer in Water. J. Chem. Phys. 2005, 122, 184506.



1

73. Resasco, D. E.; Wang, B.; Crossley, S., Zeolite-Catalysed C-C Bond Forming Reactions

2

for Biomass Conversion to Fuels and Chemicals. Catalysis Science & Technology 2016, 6,

3

2543-2559.

4

74. Wang, B.; Wright, D.; Cliffel, D.; Haglund, R.; Pantelides, S. T., Ionization-Enhanced

5

Decomposition of 2, 4, 6-Trinitrotoluene (Tnt) Molecules. J. Phys. Chem. A 2011, 115,

6

8142-8146.

7

75. Edsall, J. T.; Blanchard, M. H., The Activity Ratio of Zwitterions and Uncharged

8

Molecules in Ampholyte Solutions. The Dissociation Constants of Amino Acid Esters. J.

9

Am. Chem. Soc 1933, 55, 2337-2353.

10

76. Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J.

11

R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A., Gaussian 09 Revision A.1.

12

Gaussian Inc. 2009.

13 14

77. Head-Gordon, M.; Pople, J. A.; Frisch, M. J., MP2 Energy Evaluation by Direct Methods. Chem. Phys. Lett. 1988, 153, 503-506.

15

78. Stephens, P.; Devlin, F. J.; Chabalowski, C. F.; Frisch, M. J., Ab-Intio Calculation of

16

Vibrational Absorption and Circular-Dichroism Spectra Using Density-Function

17

Force-Fields. 1994.

18

79. Zhao, Y.; Truhlar, D. G., The M06 Suite of Density Functionals for Main Group

19

Thermochemistry, Thermochemical Kinetics, Noncovalent Interactions, Excited States,

20

and Transition Elements: Two New Functionals and Systematic Testing of Four

21

M06-Class Functionals and 12 Other Functionals. Theor. Chem. Acc. 2008, 120, 215-241.

22

80. Schlegel, H. B., Optimization of Equilibrium Geometries and Transition Structures. J.

23 24 25

Comput. Chem. 1982, 3, 214-218. 81. Fukui, K., The Path of Chemical Reactions - the IRC Approach. Acc. Chem. Res. 1981, 14, 471-476.


Page 32 of 35

Page 33 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1 2 3 4

82. Zhao, S.; Wu, J., An Efficient Method for Accurate Evaluation of the Site-Site Direct Correlation Functions of Molecular Fluids. Mol. Phys. 2011, 109, 2553-2564. 83. Zhao, S.; Liu, H.; Ramirez, R.; Borgis, D., Accurate Evaluation of the Angular-Dependent Direct Correlation Function of Water. J. Chem. Phys. 2013, 139, 034503.

5

84. Levesque, M.; Vuilleumier, R.; Borgis, D., Scalar Fundamental Measure Theory for Hard

6

Spheres in Three Dimensions: Application to Hydrophobic Solvation. J. Chem. Phys.

7

2012, 137, 034115.

8 9 10 11 12 13 14 15 16 17

85. Zhao, S.; Jin, Z.; Wu, J., New Theoretical Method for Rapid Prediction of Solvation Free Energy in Water. J. Phys. Chem. B 2011, 115, 6971-6975. 86. Yu,

Y.

X.; Wu,

J.,

Structures

of

Hard-Sphere

Fluids

from a

Modified

Fundamental-Measure Theory. J. Chem. Phys. 2002, 117, 10156-10164. 87. Carnahan, N. F.; Starling, K. E., Equation of State for Nonattracting Rigid Spheres. J. Chem. Phys. 1969, 51, 635-636. 88. Byrd, R. H.; Lu, P.; Nocedal, J.; Zhu, C., A Limited Memory Algorithm for Bound Constrained Optimization. Siam Journal on Scientific Computing 1995, 16, 1190-1208. 89. Frigo, M.; Johnson, S. G., The Design and Implementation of FFTW3. Proceedings of the IEEE 2005, 93, 216-231.

18

90. Berendsen, H. J. C.; Postma, J. P. M.; Gunsteren, W. F. V.; Hermans, J., Interaction

19

Models for Water in Relation to Protein Hydration; Springer Netherlands, 1981, p

20

331-342.

21

91. Breneman, C. M.; Wiberg, K. B., Determining Atom-Centered Monopoles from

22

Molecular Electrostatic Potentials. The Need for High Sampling Density in Formamide

23

Conformational Analysis. J. Comput. Chem. 1990, 11, 361-373.



1

92. Jorgensen, W. L.; And, D. S. M.; Tiradorives, J., Development and Testing of the OPLS

2

All-Atom Force Field on Conformational Energetics and Properties of Organic Liquids. J.

3

Am. Chem. Soc. 1996, 118, 11225-11236.

4 5

93. Barone, V.; Adamo, C.; Lelj, F., Conformational Behavior of Gaseous Glycine by a Density Functional Approach. J. Chem. Phys. 1995, 102, 364-370.

6

94. Nguyen, D. T.; Scheiner, A. C.; Andzelm, J. W.; Sirois, S.; Salahub, D. R.; Hagler, A. T., A

7

Density Functional Study of the Glycine Molecule: Comparison with Post-Hartree-Fock

8

Calculations and Experiment. J. Comput. Chem. 2015, 18, 1609-1631.

9

95. Császár, A. G., Conformers of Gaseous α-Alanine. J. Phys. Chem. 1996, 100, 3541-3551.

10

96. Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L.,

11

Comparison of Simple Potential Functions for Simulating Liquid Water. J. Chem. Phys.

12

1983, 79, 926-935.

13 14

97. Haberfield, P., What Is the Energy Difference between H2NCH2CO2H and +H3NCH2CO2-? J. Chem. Educ. 1980, 57, 346-347.

15


Page 34 of 35

Page 35 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1 2

TOC Graphic


Development of Reaction Density Functional Theory and Its

Recommend Documents