Hello,

I have been using PSI4 to investigate non-covalent interactions for associating systems (hydrogen bonding and mixed influence). All of the calculations I have done so far have been in vacuum, but most of the applications I am interested in are in solution. To get an understanding of behavior in solution, I have started to explore PCMSolver. The first trial system was the water-water complex and my work flow was as follows.

- Optimize geometry in vacuum (MP2 with aug-cc-pVTZ).
- Calculate thermodynamic properties in vacuum.
- Re-optimize geometry in solution (MP2/PCM with aug-cc-pVTZ).
- Calculate thermodynamic properties in solution.

The hardest step seems to be re-optimizing the geometry with PCM. I used the finite difference approach (dertype=‘energy’), and it took forever for the optimization to converge (59 iterations, even when starting with an optimized geometry in vacuum). Even though it took a while, the resulting geometry seems reasonable. The calculated hydrogen bond length shortened from 1.929 angstroms to 1.834 angstroms. The relative orientation between the water molecules also changed a bit (see the xyz coordinates shown below). Do these changes in geometry seem like a reasonable result for a hydrogen bonding complex? Should I expect bond lengths to shorten a bit for most complexes when re-optimizing the geometry in water?

In Vacuum (MP2 with aug-cc-pVTZ):

`O -1.478528905013 -0.073068979734 0.004135826010 H -1.843856649950 0.812352837817 -0.016261037228 H -0.521931038934 0.057289764657 -0.000318208491 O 1.406567514031 0.069106936961 -0.003826657324 H 1.755671715105 -0.386555474341 0.765611279704 H 1.752195469767 -0.420206634134 -0.753938765385`

In Solution (MP2/PCM with aug-cc-pVTZ, solvent was water):

`O -1.442781654230 0.039594940486 0.059699655839 H -1.685244466118 0.942498597239 -0.160539782306 H -0.472335854044 0.031272046541 0.016166508974 O 1.358465075784 -0.042253313013 -0.059803175369 H 1.787088761920 -0.206022777191 0.785355001878 H 1.708656845205 -0.725557565457 -0.639338798480`

I am also interested to know if the input parameters I used configure PCMSolver are reasonable for investigating hydrogen bonded geometries in solution. Most parameters have been left at the default values suggested in the PCMSolver documentation. However, I did change ‘RadiiSet = UFF’ and left ‘Scaling = False’. What is the impact of setting ‘Scaling = True’? Sorry if these are basic questions, but I am just beginning to learn about the Polarizable Continuum Model. Any references you can point me to are appreciated.

set {

basis aug-cc-pVTZ

opt_coordinates cartesian

scf_type pk

memory 16000mb

pcm true

pcm_scf_type total

}pcm = {

Units = Angstrom

Medium {

SolverType = IEFPCM

Solvent = Water

}Cavity {

RadiiSet = UFF

Type = GePol

Scaling = False

Area = 0.3

Mode = Implicit

}

}optimize(‘MP2’, molecule=water_water, dertype=‘energy’)

I know this is a long post, but I’m really trying to understand the basics before moving ahead to the other associating complexes in my investigation (30 systems in total). Thanks in advance for any assistance the community can provide. It is very much appreciated!

John