Ab Initio Electronic Structure Calculation of the ... - ACS Publications

Linda Y. Zhang and Richard A. Friesner*. Department of Chemistry, Columbia University, New York, New York 10027. Received: March 14, 1995; In Final Fo...
0 downloads 0 Views 425KB Size
J. Phys. Chem. 1995,99, 16479-16482

16479

Ab Initio Electronic Structure Calculation of the Redox Potentials of Bacteriochlorophyll and Bacteriopheophytin in Solution Linda Y. Zhang and Richard A. Friesner* Department of Chemistry, Columbia University, New York, New York 10027 Received: March 14, 1995; In Final Form: July 21, 1 9 9 P

We have carried out self-consistent reaction field (SCRF) calculations on the bacteriochlorophyll (BChl) and bacteriopheophytin (BPh) molecules at the Hartree-Fock level using the a b initio electronic structure code PSGVB together with an accurate numerical Poisson-Boltzmann solver. These are by far the largest a b initio SCRF computations carried out to date. We calculate redox potentials of these molecules in DMF solution for several geometries and compare the results with experiments of Fajer and co-workers. While the absolute redox energies are substantially in error, due to neglect of electron correlation, the relative redox energies agree qualitatively with the experimental findings. Treatment of solvation is crucial in predicting the proper ordering of the BChl and BPh redox energies. These results have implications for calculation of redox energies of the chromophores in the photosynthetic reaction center.

Introduction The functioning of the primary electron transfer in the bacterial photosynthetic reaction center is critically dependent upon the energetics of the neutral and ionic states of the bacteriochlorophyll (BChl) and bacteriopheophytin (BPh) chromophores that serve as electron donors and acceptors.I6 There have been numerous attempts to calculate the relative energies of these states, typically relying on empirical modelingI2 or semiempirical quantum chemistry.' However, considerable controversy still exists concerning the redox potentials of various species (e.g. the state P+B-, where P is the special pair primary donor and B is the intermediate BChl), the values of which are essential to understanding the charge separation mechanism. The application of ab initio quantum chemistry at a suitable level in principle is capable of resolving these issues; however, the large size of the relevant molecular structures and the necessity for considering a complex, condensed phase environment have up to now discouraged such studies. During the past several years, we have been developing an ab initio quantum chemical code, PSGVB, which is capable of handling these sorts of p r o b l e m ~ . ~ % In~particular, -~ PSGVB has efficient numerical algorithms for solving self-consistent field (SCF) equations such as the Hartree-Fock (HF),'generalized valence bond (GVB),8 and local MP29 equations for large systems, as well as an accurate self-consistent reaction field (SCRF) numerical method,I5 based upon the Poisson-Boltzmann formalism,I3 for treating solvation effects. For a series of small molecules, the SCRF methodology yielded solvation free energies with an average error of 0.6 kcallmol at the GVB level and 1.1 kcallmol at the HF leveL6 In the present paper, we present a first attempt to apply the SCRF method, at the Hartree-Fock level, to the BPh and BChl chromophores. In order to assess the accuracy of the method, we calculate redox potentials of these porphyrins in organic solvent and compare the results with experiments of Fajer and co-workers. We examine the effects of utilizing different geometries for the BChl molecules, taken from the bacterial reaction center. The calculation of electron affhities, even in the gas phase, is notoriously difficult, and the Hartree-Fock level of theory is not really adequate to determine absolute redox energies with @

Abstract published in Advance ACS Abstracts, October 1, 1995.

'CH,

Figure 1. Structure of BPh (bacteriopheophytin)used in the present work.

precision. For relative energies, however, one might hope to obtain at least qualitatively reasonable results. At a minimum, these calculations can tell us how large the solvation effects are and can allow us to make an initial assessment of the effects of geometrical distortions on the relative redox energies. The results obtained below, while not yet quantitative, are quite interesting in this regard.

Methods The ab initio electronic structure code PSGVB was utilized for all calculations, and computations were performed on an IBM 580 workstation. The calculations were carried out using the 6-31G** basis set of Pople and co-workers,I0which leads to a basis set size of 880 for BChl and 889 for BPh. Crystal structures of BChl and BPh were obtained from the Brookhaven data bank, using the X-ray coordinates from the photosynthetic reaction center structure of Diesenhofer, Michl, and coworkers.I8.l9 Initial experiments indicated that removal of the phytyl tail did not substantially affect the charge distributions, so the calculations below are actually on structures 1 and 2, displayed in Figures 1 and 2, respectively. The experimental results for redox potentials have been reported in the solvent DMF.'a2 To model this solvent in our SCRF approach, we must specify its dielectric constant and the size of a solvent molecule (approximated as spherical), which is used to generate the Connally surface for the dielectric interface.I3 The dielectric constant of DMF is 36.7, and the radius of the probe solvent molecule is taken to be 2.67 A.

0022-365419512099-16479$09.00/0 0 1995 American Chemical Society

Zhang and Friesner

16480 J. Phys. Chem., Vol. 99, No. 44, 1995

and the values of (qglH0lqg) (the gas phase energy), { ( ~ s l H ~ l q s ) I/2(qsIHept Hnpll~s)} (the free energy in solution), and AGs (the solvation energy) were presented in Table 1. The key parameters in this SCRF approach are the dielectric radii of the atoms, which define the dielectric interface via the Connally surface. We have previously parametrized radii for C, 0, N, and H,6 and we did not readjust these values in the present work for both neutral and ionic species. The dielectric radius of the Mg atom in BChl was set to 1.51, which was taken from the DREIDING force field of Goddard and co-workers.20 We note that since the solvent is different (previous work investigated aqueous solvation) and neutral and ionic species are different (previous work only investigated neutral molecules), one might expect that reparametrization of the radii would lead to quantitatively improved results. This question will be investigated further in future publications. The full solvation free energy should include the solute cavity energy, which is typically estimated as being proportional to molecular surface area. However, only the difference of the solvation energies of the neutral and ionic molecules contributes to the redox energy, which is what we are interested in the present work. Because the same atomic radii were used for neutral and ionic species, the cavity energies are equal; thus, the cavity energies do not affect the redox energies. Absolute Redox Energy Calculations. The reduction of a molecule M in solution can be expressed as

+

+

Figure 2. Structure of BChl (bacteriochlorophyll) used in the present work. BChl-I and BChl-SP have the same structures but different threedimensional geometries; BChl-I is the intermediate BChl on the L-side of the bacterial photosynthetic reaction center, and BChl-SP is one BChl of the special pair.

SCRF Calculations. The details of our SCRF methodology have been described previously in ref 6. Considering the solvent as a dielectric continuum, we employ the program DelPhi to solve the relevant Poisson e q ~ a t i 0 n . l To ~ briefly review the approach, the gas phase wave function is initially generated in PSGVB, fit to atomic point charges via electrostatic potential (ESP) fitting, and then passed to DelPhi to generate a reaction field, which is represented by another set of point charges distributed on the surface of the solute molecule. The quantum chemistry is now carried out in the presence of the reaction field, leading to polarization of the solute and determination of a new set of point charges. This process is repeated until the solvation energy and quantum chemical wave function have converged. The solvation energy here (as presented in Table 1) is the electrostatic free energy of solvation. The Hamiltonian of the system in the reaction field could be written as (Ho Hept Hnpt),where HO is the gas phase Hamiltonian of the solute molecule, Heptis the sum of Coulomb interactions between electrons and those surface point charges (the reaction field), and Hnptis the sum of Coulomb interactions between the solute nuclear charges and those surface point charges. Then the total quantum mechanical energy of the solute phase is given by (qSIHo Hept Hnptlqs), where qSis the solvated solute wave function. But those surface charges are virtual, representing the reaction field of the fully relaxed dielectric continuum; consequently the (free) energy associated with this polarization should be 1/2(qSIHept Hnptlqs)(see refs 6, 22, and 23 for detailed discussions). So the final expression of the solvation energy is given by

+

+

-

M-(s)

AGS(M-M-)

(3)

where AGS(M-M-) denotes the absolute redox energy of the above process. The free energy of the electron at rest in the gas phase is set to zero.21 Then, the redox energy can be calculated via a thermodynamic cycle, where AGS(M) and

+

+

+

M(s) -te-

1

M(g)

AB(M - M )

1

Aff(W M-(g)

AGs(M-) are the solvation energies of M and M-, respectively, and hEg(M-+M-) is the gas phase energies difference between M and M-, defined as the gas phase redox energy . We can obtain the absolute redox energy AGs(M-M-) as AGS(M-M-)

+ AGS(M) - AGs(M-)

= A,??(M-M-)

(4)

So by calculating the gas phase energies and solvation energies of M and M-, respectively, we can get the absolute redox energy of M in solution. Results

where qsis the gas phase wave function of the solute molecule

Table 1 presents gas phase energies, free energies in solution, and solvation energies for BChl and BPh in both the neutral and anionic forms for several molecular geometries. In particular, we have carried out calculations for a BChl geometry of a molecule in the special pair (designated BChl-SP in Table 1) and for the geometry of the intermediate bacteriochlorophyll

TABLE 1: Calculated Gas Phase Energies, Free Energies in Solution, and Solvation Energies species BPh BChl-I BChl-SP a

gas phase energy (hartrees) neutral anionic neutral anionic neutral anionic

1 hartree = 627.5059 kcal/mol.

-2052.203 -2052.250 -2250.839 -2250.891 -2250.864 -2250.915

179 727 98 780 5 19 36 941 257 68 940 699 97 797 614 76 577 301 83

free energy in DMF (hartrees) -2052.256 346 543 90 -2052.353 900 644 53 -2250.913 145 260 89 -225 1.008 498 305 95 -2250.931 347 015 08 -225 1.029 761 406 60

solvation energy in DMF (kcaVmo1") -33.362 -64.708 -45.936 -73.141 -41.760 -71.651

682 074 34 858 183 58 207 455 08 005 050 14 380 923 32 61049656

Redox Potentials of Bacteriochlorophyll and -pheophytin

J. Phys. Chem., Vol. 99,No. 44, 1995 16481

TABLE 2: Calculated Redox Energies and Experimental Results PSGVB (6-31G**) experiment" gas phase in DMF in DMF (kcaVmo1) (kcaVmo1) (eV) (kcaVmo1) (eV) BPh 29.8698 61.2160 2.655 92.242 4.000 BChl-I 32.6303 59.8351 2.595 87.630 3.800 BChl-SP 31.8645 61.7557 2.678 diff(BPh - BChl-I) -2.7605 1.3809 0.060 4.612 0.200 These experimental data are absolute reduction potentials converted from those vs NHE (normal hydrogen electrode) by eq 7. (I

(designated BChl-I in Table 1). From these results , we compute the gas phase and solution phase redox energies of all three molecules (Table 2). The following qualitative features are apparent: (1) Solvation makes a major contribution to the absolute redox potential of both species, on the order of 1 eV. (2) The difference in redox potential between BPh and BChl has contributions from both the gas phase and solvation terms. We next compare our results for the solution phase redox potentials with the experiments of Fajer and co-workers. The latter are reported vs the NHE (normal hydrogen electrode) and hence have to be converted to absolute redox potentials to make comparisons with the quantum chemical calculations. NHE, which is also called the absolute potential of the hydrogen electrode, is the free energy change E", of the process

The best estimate for ,$!, is 4.5 eV.,' Reactions 3 and 5 can be coupled as

L? = AG'(M-M-)

- J$ = AG'(M-M-)

- 4.50 eV

(7)

is the relative reduction potential for M, which can be measured by experiment. By using eq 7, we converted the raw experimental data, 500 mV for BPh and 700 mV for BChl, into the absolute values 4.00 V and 3.80 V, respectively, and compared them to our calculations in Table 2. We also compute the difference in redox energy for each of the BChl geometries as compared to the BPh molecule, in both the gas phase and solution phase. These results are presented in Table 2 as well. The first conclusion from the above comparison is that the Hartree-Fock level of theory is inadequate to calculate electron affinities; the errors in the present case are on the order of 1 eV. This is unsurprising because the ionic species will have a much larger correlation energy than the neutrals, due to the presence of the extra electron. The HF level calculations duly underestimate the electron affinity, and the magnitude of the error is consistent with the results of other such calculations on other molecules. Despite this sizable correlation error, it is possible that relative redox potentials might be rendered more accurately. The first interesting observation concerning the differences in redox potential is that there are significant effects of altering the geometry. In fact, there is a sign change as one goes from BChlSP (which is more easily reduced than BPh) to BChl-I (which is more difficult to reduce than BPh). It is possible that these effects contribute to the function of BChl-SP in the special pair, where the internal charge-separated state of the dimer has been shown by our group, from an analysis of spectroscopic data, to lie at a very low en erg^;^.^ however, verifying this hypothesis

will require investigation of the SP in the reaction center environment. We expect that the BChl-SP geometry is substantially more distorted from the geometry of a molecule in DMF than that of BChl-I, as the structure of the former is undoubtedly influenced strongly by the close-lying second molecule of the dimer. Thus, the BChl-I results should be closest to the experimental values obtained in DMF. This is indeed the case; BChl-I is predicted to be 1.8 kcal/mol harder to reduce than BPh while the experimental value is -4.6 kcdmol. This agreement is not quantitative, but it is encouraging in that the correct sign and order of magnitude were obtained. An important observation is that the solvation term actually reverses the order of redox energies for the two molecules, thus emphasizing that, to understand the energetics of BChl molecules in the reaction center, accurate modeling of the environment is essential.

Conclusion We have shown in this paper that ab initio SCRF methods at the Hartree-Fock level are capable of yielding reasonable qualitative predictions of the relative redox potentials in solution for the BChl and BPh chromophores. In view of the approximations inherent in the calculation and of the uncertainty of the experiments, the level of agreement obtained is satisfactory for the present. When electron correlation is added to the calculations (we intend to carry out local MP2 calculations as the next level of theory9 and follow these up with GVB-RCIMP2 calculations which should be quantitatively accurate) and parametrization of the SCRF model becomes more accurate, one would want to do significantly better; it remains to be seen whether this will be the case. Our intention, after carrying out correlated calculations on BChl and BPh in solution, is to apply this methodology to the photosynthetic reaction center, treating most of the protein as a dielectric continuum in which point charges are embedded, and to nearest neighbor residues at the explicit quantum chemical level. If the use of correlation yields the expected increase in accuracy, real progress can be made in evaluating the energy of the intermediate BChl charge separated state, a key quantity in understanding the mechanism of charge separation. The variation in relative energies as a function of geometry and environment in the present paper should serve as a cautionary note in the evaluation of many of the previous electronic structure computations on the reaction center. In general these methods have been semiempirical in character and the treatment of solvation has been either nonexistent or at a relatively crude level with regard, for example, to treatment of the molecular cavity shape. In no cases were accurate quantum chemical geometry optimizations in the reaction center carried out. In view of our results, it is thus unlikely that these calculations had an accuracy better than -0.2 eV (of course, the accuracy could have been worse than this for unrelated reasons, such as inaccuracies in the semiempirical parametrization). There is obviously much to be learned from semiempirical studies even if this level of accuracy was not achieved; however, some of the questions in the reaction center, such as the relative energies of the excited and charge-separated state of the special pair, will require this sort of precision. It will take a lot of work to produce results of this quality for the reaction center, and it is our intention to proceed along this path despite the obvious difficulties.

Acknowledgment. This work was supported in part by grants to R.A.F. from the NIH and the NSF. We thank Marshall Newton for a number of useful discussions.

Zhang and Friesner

16482 J. Phys. Chem., Vol. 99, No. 44, 1995 References and Notes (1) Fajer, J.; Davis, M. S.; Brune, D. C.; Spaulding, L.; Borg, D.; Forman. A. Brookhaven Svma Biol. 1976, 28. 74, (2) Davis, M. S.; Forkan, A.; Hanson. L. K.; Fajer, J. J. Phys. Chem. 1979, 83, 3325. ( 3 ) Friesner, R. Annu. Rev. Phys. Chem. 1991, 42, 341. (4) Lathrop, E.; Friesner, R. J. Phys. Chem. 1994, 98, 3056. ( 5 ) Lathrop, E.; Friesner, R. J. Phys. Chem. 1994, 98, 3050. (6) Tannor, D.; Marten, B.; Murphy, R.; Ringnalda, M.; Friesner, R.; Nicholls, A.; Goddard, W. A., 111; Honig, B. J. Am. Chem. SOC. 1994, 116, 11875. (7) Greeley, B.; Russo, T.; Mainz, D.; Friesner, R.; Langlois, J. M.; Goddard,W., 111; Ringnalda, M. J. Chem. Phys. 1994, 101, 4028. (8) Murphy, R.; Friesner, R.; Goddard, W., 111; Ringnalda, M. J. Chem. Phys. 1994, 101, 2986. (9) Murphy, R.; Beachy, M.; Ringnalda, M.; Friesner, R. J. Chem. Phys., submitted. (IO) Hehre, W.; Ditchfield, R.; Pople, J. A. J. Chem. Phys. 1972, 56, 2257. Harihan, P. C.; Pople, J. A. Theor. Chim. Acta (Berlin) 1973, 28, 213. (11) Thompson, M. A.; Zemer, M. C. J. Am. Chem. SOC. 1991, 113, 8210.

(12) Bixon, M.; Jortner, J.; Michl-Beyerle, M. E. Biochim. Biophys. Acta 1991, 1056, 301.

(13) Honig, B.; Sharp, K.; Yang, A. J . Phys. Chem. 1993, 97, 1101. (14) Nicholls A,; Honig, B. J. Comput. Chem. 1991, 12, 435. (15) Tomasi, J.; Algona, R.; Bonnacorsi, R.; Ghio, C. In Modeling of Structures and Properties of Molecules; Maksic, C., Ed.; Ellis Horwood: Chichester, U.K., 1987; p 330. (16) Friesner, R.; Won, Y. Biochim. Biophy. Acta 1989, 977, 99. (17) Hansen, L. K. Photochem. Photobiol. 1988, 47, 903. (18) Deisenhofer, J.; Michl, H. Science 1989, 245, 463. (19) Deisenhofer, J.; Epp, 0.;Huber, R.; Michl, H. Nature 1985, 318, 618. (20) Mayo, S.; Olafson, B.; Goddard, W. A,, 111 J. Phys. Chem. 1990, 94, 897. (21) Pearson, It. J. Am. Chem. Soc. 1986, 108, 6109. (22) Jackson, J. D. Challical Electrodynamics, 2nd ed.; Wiley: New York, 1975; pp 158-162. (23) Newton, M. D. J. Chem. Phys. 1973, 58, 5833. Newton, M. D. J. Phys. Chem. 1975, 79, 2795.

JP9507253