<jupyter_start><jupyter_text>Tutorial 11: Intermolecular interaction energies (SAPT) and solvation effects (PCM)In this tutorial you will explore:- how to compute molecular interation energies with the SAPT0 method- how. to account for solvation effects via the polarizable continuum model (PCM)<jupyter_code>import psi4<jupyter_output><empty_output><jupyter_text>1. Analysis of the intermolecular interaction in the HF dimer with SAPT0In this first example, we will consider the interaction energy of the HF dimer. To compute the interaction energy, we optimize the geometry of the dimer with DFT and then do a SAPT0 computation. Geometry optimizationFor weakly bound dimers like (HF)$_2$ it is important to optimize the geometry using both Cartesian and internal coordinates, and to compute the full hessian at the beginning to speed up the convergence of the geometry optimization.The optimization was performed with the following input```pythonhf_dimer = psi4.geometry("""0 1F       -1.3976919319      0.1218846771      0.0000000000                 H       -0.8539861467     -0.7000656533      0.0000000000                 F       -0.3234710479     -2.9923356673      0.6218021310                 H       -0.2279209870     -3.9731967777      0.6218021310                  units angstrom no_reorient symmetry c1""")psi4.set_output_file('hf_dimer_opt.txt', False)psi4.set_options({'OPT_COORDINATES' : 'BOTH', 'FULL_HESS_EVERY' : 0})psi4.optimize('b3lyp/aug-cc-pvdz')```and it yielded the following optimized geometry```pythonF                    -0.597917227851     1.185245761736    -0.435087735701H                    -0.132328360109     0.460602356132    -0.075129151973F                     0.592935186451    -1.111803010103     0.438932749487H                     0.226244295964    -1.845063904792     0.002647203482```  SAPT0 ComputationWe can now assemble the input for the SAPT0 computation. To run SAPT, we need to specify the geometries of the two interacting molecules.  In Psi4, you can use the `--` separator to tell psi4 which atoms belong to the first and second molecule in the `psi4.geometry` function. For each fragment you can also specify the charge and the multiplicity. Note that SAPT does not use spatial symmetry (so it is a good idea to specify `symmetry c1` in the input) and **needs the geometry of the system to remain fixed throughout monomer and dimer calculations** (by adding the `no_reorient` keyword at the end of the input).<jupyter_code># we take the optimized geometry and use the -- separator to tell Psi4 that
# the first two atoms belong to molecule A and the last two to molecule B
hf_dimer = psi4.geometry(
"""
0 1
F                    -0.597917227851     1.185245761736    -0.435087735701
H                    -0.132328360109     0.460602356132    -0.075129151973
--
0 1
F                     0.592935186451    -1.111803010103     0.438932749487
H                     0.226244295964    -1.845063904792     0.002647203482              
 no_reorient
 symmetry c1
"""
)

# If you have fortecubeview installed, you can visualize the geometry using the following two lines
# import fortecubeview
# fortecubeview.geom(molecule=hf_dimer)<jupyter_output><empty_output><jupyter_text>Next we can run a SAPT0 computation. The input follows the standard format `psi4.energy('sapt0/')` where `` is the label of the basis set (e.g., def2-SVP). This function returns the total SAPT0 **interaction energy** in units of Hartree. A negative number indicates that the complex of the two molecules is more stable than the individual separated molecules.<jupyter_code>psi4.set_output_file('hf_dimer_sapt0.txt', False)
psi4.energy('sapt0/aug-cc-pvdz')<jupyter_output><empty_output><jupyter_text>To see the energy broken down into all of its component you can open the output file or print the content in a cell with the command:```!cat hf_dimer_sapt0.txt```Psi4 prints a table in the output file that shows all the contributions to the SAPT0 energy:To understand the meaning of all of these corrections, I recommend that you take a look at the "Psi4 SAPT manual page (http://www.psicode.org/psi4manual/master/sapt.html) ```    SAPT Results   --------------------------------------------------------------------------------------------------------    Electrostatics                -11.30028263 [mEh]      -7.09103440 [kcal/mol]     -29.66888795 [kJ/mol]      Elst10,r                    -11.30028263 [mEh]      -7.09103440 [kcal/mol]     -29.66888795 [kJ/mol]    Exchange                        8.93916328 [mEh]       5.60940964 [kcal/mol]      23.46976995 [kJ/mol]      Exch10                        8.93916328 [mEh]       5.60940964 [kcal/mol]      23.46976995 [kJ/mol]      Exch10(S^2)                   8.89619835 [mEh]       5.58244875 [kcal/mol]      23.35696556 [kJ/mol]    Induction                      -3.71576506 [mEh]      -2.33167778 [kcal/mol]      -9.75573983 [kJ/mol]      Ind20,r                      -4.57578012 [mEh]      -2.87134537 [kcal/mol]     -12.01370904 [kJ/mol]      Exch-Ind20,r                  2.13148785 [mEh]       1.33752882 [kcal/mol]       5.59622058 [kJ/mol]      delta HF,r (2)               -1.27147280 [mEh]      -0.79786123 [kcal/mol]      -3.33825137 [kJ/mol]    Dispersion                     -2.10841994 [mEh]      -1.32305349 [kcal/mol]      -5.53565579 [kJ/mol]      Disp20                       -2.49992880 [mEh]      -1.56872901 [kcal/mol]      -6.56356216 [kJ/mol]      Exch-Disp20                   0.39150886 [mEh]       0.24567552 [kcal/mol]       1.02790637 [kJ/mol]      Disp20 (SS)                  -1.24996440 [mEh]      -0.78436450 [kcal/mol]      -3.28178108 [kJ/mol]      Disp20 (OS)                  -1.24996440 [mEh]      -0.78436450 [kcal/mol]      -3.28178108 [kJ/mol]      Exch-Disp20 (SS)              0.23971773 [mEh]       0.15042515 [kcal/mol]       0.62937881 [kJ/mol]      Exch-Disp20 (OS)              0.15179113 [mEh]       0.09525037 [kcal/mol]       0.39852756 [kJ/mol]  Total HF                         -6.07688441 [mEh]      -3.81330254 [kcal/mol]     -15.95485783 [kJ/mol]  Total SAPT0                      -8.18530435 [mEh]      -5.13635603 [kcal/mol]     -21.49051362 [kJ/mol]  Special recipe for scaled SAPT0 (see Manual):    Electrostatics sSAPT0         -11.30028263 [mEh]      -7.09103440 [kcal/mol]     -29.66888795 [kJ/mol]    Exchange sSAPT0                 8.93916328 [mEh]       5.60940964 [kcal/mol]      23.46976995 [kJ/mol]    Induction sSAPT0               -3.68473308 [mEh]      -2.31220492 [kcal/mol]      -9.67426537 [kJ/mol]    Dispersion sSAPT0              -2.10272003 [mEh]      -1.31947674 [kcal/mol]      -5.52069067 [kJ/mol]  Total sSAPT0                     -8.14857246 [mEh]      -5.11330642 [kcal/mol]     -21.39407405 [kJ/mol]  --------------------------------------------------------------------------------------------------------```  2. Modeling the solvation energy of ethanol in waterIn this example, we are going to use the PCM model to evaluate the electrostatic contributions to the solvation energy of ethanol. We use a geometry optimized at the B3LYP/def2-SVP level of theory in the absence of water molecules (gas phase).First, we compute the energy of ethanol in the gas phase. Note that we added the strings```no_reorientno_comsymmetry c1```These prevent Psi4 from reorienting the geometry, shifting it to the center of mass, and enforce the use of no symmetry. In this case these keywords do not make any difference, but it is good practice to keep them here so thatn when we switch on PCM everything will be consistent.<jupyter_code>ethanol_geom = psi4.geometry("""
    C           -0.333009979899    -0.423615357612    -0.184200903142
    C            1.179153362301    -0.414706417524    -0.313836658083
    H           -0.651311352848    -1.316907925328     0.393604041781
    H           -0.791894825341    -0.517368555358    -1.191066531406
    O           -0.732927481926     0.773239648281     0.455152339143
    H            1.537340993654    -1.331700551471    -0.808276264206
    H            1.649040960998    -0.346302943607     0.679936874753
    H            1.508378518884     0.453776056275    -0.905787283729
    H           -1.694347581166     0.768383333798     0.538039164456
no_reorient
no_com
symmetry c1
""")

psi4.set_output_file('ethanol_gas.txt', False)
ethanol_gas = psi4.energy('b3lyp/def2-SVP')
print(f'Energy of ethanol in the gas phase: {ethanol_gas}')<jupyter_output><empty_output><jupyter_text>In the next step we setup a PCM computation. The most important option here is the `Solvent` keyword, which is set to `Water`. Other options are:```Water , H2O;Propylene Carbonate , C4H6O3;Dimethylsulfoxide , DMSO;Nitromethane , CH3NO2;Acetonitrile , CH3CN;Methanol , CH3OH;Ethanol , CH3CH2OH;Acetone , C2H6CO;1,2-Dichloroethane , C2H4CL2;Methylenechloride , CH2CL2;Tetrahydrofurane , THF;Aniline , C6H5NH2;Chlorobenzene , C6H5CL;Chloroform , CHCL3;Toluene , C6H5CH3;1,4-Dioxane , C4H8O2;Benzene , C6H6;Carbon Tetrachloride , CCL4;Cyclohexane , C6H12;N-heptane , C7H16;```The other options just repeat default values in the PCM module and you should not tweak them. To pass the PCM options you can call the function `psi4.pcm_helper(pcm_string)` passing the `pcm_string` that we define on top.Please consult the PCM manual page (https://pcmsolver.readthedocs.io/en/latest/users/input.htmlmedium-section-keywords) for more information on how to run PCM computations.<jupyter_code>pcm_string = """
    Units = Angstrom
    Medium {
        SolverType = IEFPCM
        Solvent = Water
    }
    Cavity {
       RadiiSet = UFF
       Type = GePol
       Scaling = False
       Area = 0.3
       Mode = Implicit
    }
"""

# set the PCM options
psi4.pcm_helper(pcm_string)

# tell Psi4 to run a PCM computation
psi4.set_options({'pcm': True, 'pcm_scf_type': 'total'})

psi4.set_output_file('ethanol_aq.txt', False)
ethanol_aq = psi4.energy('b3lyp/def2-SVP')

print(f'Energy of ethanol in water (PCM):   {ethanol_aq}')
print(f'Solvation energy: {627.51 * (ethanol_aq - ethanol_gas):.2f} kcal/mol')<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.10653057179665049, 'fraction_numerical': 0.1559805767515608, 'mean_word_length': 2.4060074249071888, 'pattern_counts': {'":': 0, '<': 21, '<?xml version=': 0, '>': 21, 'https://': 1, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 2, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '5', 'dir': 'jupyter-structured-clean-dedup', 'id': '458022', 'license': 'BSD-3-Clause', 'n_tokens_mistral': 4381, 'n_tokens_neox': 3366, 'n_words': 970, 'path': 'Example/CHEM371/Notebooks/11-ModelingInteractions/11-ModelingInteractions.ipynb', 'repo_name': 'yychuang/109-2-compchem-lite', 'size': '10091'}