excited-state solvent models and a toy system

My latest projects are focused around a solvent model I extended and interfaced to the excited state methods developed in our group during my PhD. But before I write about the toy system I came up with to test my latest work, let me briefly introduce you to polarizable-continuum solvation models. PCMs basically compute a set of point charges on the molecular surface which depend on the molecular electrostatic potential (ESP) to mimic the interaction with solvent molecules. These (point) charges are added to the one-electron Hamiltonian (just like the nuclei) during the self-consistent field (SCF) calculation and updated in every iteration. Eventually, one obtains a wavefunction of a molecule in solution and the interaction between the point charges and the molecular ESP corresponds to the solute-solvent interaction. The single parameter used to describe the solvent (more parameters are used, e.g. for the cavity construction) is its macroscopic dielectric constant epsilon.

Visualization of the PCM surface-charges (blue -, red +) obtained for my beloved example nitrobenzene in its electronic ground state. As one would expect, the surface-charge close to the negatively charged nitro group (front right) is positive - just like the potential nitrobenzene would exhibit in a polar solvent.

 After finishing and extensively evaluating a so-called non-equilibrium (NEq) solvent model for vertical excitation energies against a set of experimental data for nitroaromatics (find it here), I have recently been working on the respective equilibrium variant. In a nutshell, non-equilibrium and equilibrium solvent models can be rationalized in analogy to the Franck-Condon (FC) principle:

• The non-equilibrium case applies to vertical (i.e. very fast) processes like absorption and fluorescence. To be in accordance with the FC principle, you want to clamp the solvents nuclei, while relaxing its electronic degrees of freedom. Within a PCM, this is realized by using the optical dielectric constant (n²) instead of epsilon for a recomputation of the point charges w.r.t. the excited-state ESP. n² is nothing but the macroscopic polarizability ofa solvent at the frequency of visible light, whereas epsilon is the total polarizability (including a rearrangement of the solvent nuclei).
• The equilibrium case applies to long-lived states, i.e. whenever you expect your systems nuclei to relax w.r.t. an excited state potential energy surface. Consequently, you should also relax the solvents nuclei, which is equivalent to recompute the surface charge for the new molecular charge distribution using the full dielectric constant just like for the ground state.

Schematic visualization of equilibrium and non-equilibirum solvation for the example of an excited-state proton transfer in quinolines. The proton transfer can be regarded as the molecules geometrical relaxation, which is accompanied by the solvent equilibration.

Although the explanation and also the equations for the non-equilibrium case are slightly more complicated due to the separation into fast (n²) and slow (epsilon) parts of the polarization, they are much simpler to implement. The reason being that n² of typical solvents is rather small one hence, one can apply perturbation theory an turn it into a correction for excitation energies. For details, see my paper.

The equilibrium case looks straightforward at first glance, but turns out much more complicated w.r.t. its implementation because of the way our excited-state methods work: We compute the Hartree-Fock (HF) ground-state wavefunction (in the form of molecular orbitals, short MOs) and subsequently use it as a basis for the calculation of excited-state wavefunctions and excitation energies with ADC (Algebraic-Diagrammatic Construction for the Polarization Propagator, a configuration-interaction type excited-state method).
Since the HF MOs contain the interaction with PCM surface charges just like they contain the interaction with the nuclei, one needs the PCM surface charges of the excited state before you do the initial HF calculation. Since this is not possible, one has to iterate between HF and ADC until the solvent field and the excited-state ESP are converged, like so:

1. Standard HF calculation with PCM solvent model => MOs
2. ADC calculation using HF MOs => excited-state density
3. HF calculation in the "frozen" surface-charges for ES density=> new MOs
4. ADC calculation with new MOs to obtain new excited-state density
5. repeat from 3 until surface charges/energy are converged 

Eventually, I got the implementation working and the first numbers were looking promising. To check if the code does indeed work correctly, I came up with a simple toy system in which the solvated ground state is identical to an equilibrium solvated charge-transfer excited state: an asymmetrically charged, cationic ethene dimer with a distance of 1 nm:

The cationic ethene dimer as test system. The influence from the PCM is indicated by the boldness of the cavity. As soon as the iterations of the ES reaction field are converged, the system is in the same situation as in the beginning. Also, the non equilibrium situations with charge-tranfer (CT) state the in the field of the ground state (abbreviated  "0" in the scheme) (Nr. 2) and the ground state in the field of the CT state (converged Nr. 3) should be very comparable.

There is only a one practical problem: During the solvent-equilibration iterations (described above), the charge in the ground state jumps back and forth, causing the solvent-equilibration procedure to oscillate instead of converging. Obviously, the HF-SCF calculation in the frozen surface charges (point 3) of the charge-tranfer state converges to the energetically much more favorable solution with the positive charge sitting the ethene in the negatively charged cavity (No surprise, this is exactly what the SCF algorithm is supposed to do). Employing the MOs of the previous iteration step as starting point in the subsequent SCF cycles is not enough do do the trick, the other solution seems too low energetically.

How I was nevertheless able to converge the calculations using a number of tricks will be presented along with some more details on the system and preliminary results in my next post very soon.

So long!
Jan