Calculating RESP charges without a server

This tutorial demonstrates the two-part process of calculating charges without a QCFractal server, e.g. on a shared cluster. It presumes that you are already familiar with the basic steps of using PsiRESP, such as creating a molecule or the ways that a job can be customized. If not, please see the related notebooks:

Note: you will need to have Psi4 installed in the environment that you used open the notebook, to run the jobs. If the below command fails, you either need to install Psi4, or switch to the correct environment and re-open the notebook.

If you want to calculate charges without Psi4 or a QCFractal (i.e. you have already pre-computed grids and electrostatic potentials), please see Calculating RESP charges with pre-computed ESPs.

[1]:

%%bash

psi4 --version

1.5

[2]:

import psiresp

Creating the molecule

We choose molecule with multiple potential conformers, and geometry optimization is turned on to demonstrate both the geometry optimization and ESP computation steps.

[3]:

butanol = psiresp.Molecule.from_smiles(
    "CCCCO",
    optimize_geometry=True,
    conformer_generation_options=dict(n_max_conformers=3,
                                      keep_original_conformer=False),
)

Running the job

As previously, we create a typical psiresp.Job.

[4]:

job = psiresp.Job(molecules=[butanol],
                  working_directory="no_server_working_directory")

Below we will remove any files from previous runs of this example for demonstration purposes.

[5]:

%%bash

rm -rf no_server_working_directory

Geometry optimization

Firstly, we simply call job.run(). This is expected to write Psi4 execution files for geometry optimization, and then exit Python with an error message that tells us of a script we can use to run Psi4. Note that if you are running this notebook in Jupyter, the Python exit will not occur but will only print a warning message about how to properly exit IPython, which is quite nice for demonstration purposes.

[6]:

job.run()

generate-conformers: 100%|████████████████████████| 1/1 [00:06<00:00,  6.22s/it]

An exception has occurred, use %tb to see the full traceback.

SystemExit: Exiting to allow running QM computations; commands are in no_server_working_directory/optimization/run_optimization.sh

/Users/lily/anaconda3/envs/tmpenv3.9/lib/python3.9/site-packages/IPython/core/interactiveshell.py:3377: UserWarning: To exit: use 'exit', 'quit', or Ctrl-D.
  warn("To exit: use 'exit', 'quit', or Ctrl-D.", stacklevel=1)

If we have a look at what is in the directory, an optimization subdirectory has been created, and three Psi4 job files (ending in .msgpack) have been generated. These are named with the molecular name or formula, followed by the hash of the conformer coordinates, followed by the hash of the optimization options. This allows the Psi4 files to persist between jobs, and be re-used in others.

[7]:

%%bash

tree no_server_working_directory

no_server_working_directory
└── optimization
    ├── C4H10O_221cf2ccc9bf1573c95ba89dfa5458cf47a49339_c9ce731306cb83b137c5cfd5f69a120483b61005.msgpack
    ├── C4H10O_683183ed77833f6861a934f55156ef4013feb14b_c9ce731306cb83b137c5cfd5f69a120483b61005.msgpack
    ├── C4H10O_d53a466202b42f8f92a53755db4017420d14aa65_c9ce731306cb83b137c5cfd5f69a120483b61005.msgpack
    └── run_optimization.sh

1 directory, 4 files

As directed by the SystemExit error, we can run the run_optimization.sh file. The commands are quite simple:

[8]:

%%bash

cat no_server_working_directory/optimization/run_optimization.sh

#!/usr/bin/env bash
psi4 --qcschema C4H10O_683183ed77833f6861a934f55156ef4013feb14b_c9ce731306cb83b137c5cfd5f69a120483b61005.msgpack
psi4 --qcschema C4H10O_221cf2ccc9bf1573c95ba89dfa5458cf47a49339_c9ce731306cb83b137c5cfd5f69a120483b61005.msgpack
psi4 --qcschema C4H10O_d53a466202b42f8f92a53755db4017420d14aa65_c9ce731306cb83b137c5cfd5f69a120483b61005.msgpack

[9]:

%%bash

cd no_server_working_directory/optimization && bash run_optimization.sh

Psi4 actually writes the output back into the same .msgpack file, so not much has changed if we look at the directory structure.

[10]:

%%bash

tree no_server_working_directory

no_server_working_directory
└── optimization
    ├── C4H10O_221cf2ccc9bf1573c95ba89dfa5458cf47a49339_c9ce731306cb83b137c5cfd5f69a120483b61005.msgpack
    ├── C4H10O_683183ed77833f6861a934f55156ef4013feb14b_c9ce731306cb83b137c5cfd5f69a120483b61005.msgpack
    ├── C4H10O_d53a466202b42f8f92a53755db4017420d14aa65_c9ce731306cb83b137c5cfd5f69a120483b61005.msgpack
    ├── run_optimization.sh
    └── timer.dat

1 directory, 5 files

Calculating ESP

However, if we create a new job and run it again, PsiRESP checks the .msgpack files and finds the optimized coordinates. That allows the Job to progress to the single point stage for calculating the electrostatic potential.

[11]:

job = psiresp.Job(molecules=[butanol],
                  working_directory="no_server_working_directory")
job.run()

generate-conformers: 100%|██████████████████████| 1/1 [00:00<00:00, 9279.43it/s]

An exception has occurred, use %tb to see the full traceback.

SystemExit: Exiting to allow running QM computations; commands are in no_server_working_directory/single_point/run_single_point.sh

/Users/lily/anaconda3/envs/tmpenv3.9/lib/python3.9/site-packages/IPython/core/interactiveshell.py:3377: UserWarning: To exit: use 'exit', 'quit', or Ctrl-D.
  warn("To exit: use 'exit', 'quit', or Ctrl-D.", stacklevel=1)

This is in a new subdirectory.

[12]:

%%bash

tree no_server_working_directory

no_server_working_directory
├── optimization
│   ├── C4H10O_221cf2ccc9bf1573c95ba89dfa5458cf47a49339_c9ce731306cb83b137c5cfd5f69a120483b61005.msgpack
│   ├── C4H10O_683183ed77833f6861a934f55156ef4013feb14b_c9ce731306cb83b137c5cfd5f69a120483b61005.msgpack
│   ├── C4H10O_d53a466202b42f8f92a53755db4017420d14aa65_c9ce731306cb83b137c5cfd5f69a120483b61005.msgpack
│   ├── run_optimization.sh
│   └── timer.dat
└── single_point
    ├── C4H10O_1c5e35d36810e33e916c6644f1838761efd45846_e746222796fc2c4c5a1f896fa1cc1cefffe7044c.msgpack
    ├── C4H10O_6b12cd9fa9bd56c0692fcd6384bbc05af96c1c1e_e746222796fc2c4c5a1f896fa1cc1cefffe7044c.msgpack
    ├── C4H10O_8f405ba6ddee36ab6dcc30af38dfc4fea14c5a86_e746222796fc2c4c5a1f896fa1cc1cefffe7044c.msgpack
    └── run_single_point.sh

2 directories, 9 files

PsiRESP can also respond to jobs that have been partially completed. To demonstrate this, we only run two of the three required single points.

[13]:

%%bash

cat no_server_working_directory/single_point/run_single_point.sh

#!/usr/bin/env bash
psi4 --qcschema C4H10O_1c5e35d36810e33e916c6644f1838761efd45846_e746222796fc2c4c5a1f896fa1cc1cefffe7044c.msgpack
psi4 --qcschema C4H10O_6b12cd9fa9bd56c0692fcd6384bbc05af96c1c1e_e746222796fc2c4c5a1f896fa1cc1cefffe7044c.msgpack
psi4 --qcschema C4H10O_8f405ba6ddee36ab6dcc30af38dfc4fea14c5a86_e746222796fc2c4c5a1f896fa1cc1cefffe7044c.msgpack

[14]:

%%bash

cd no_server_working_directory/single_point

# this loop is in case the molecule optimizes to different coordinates
# on different runs, which will alter the filename hash accordingly
for msgpack in $( ls *.msgpack | tail -n 2 ) ; do
    echo "executing $msgpack"
    psi4 --qcschema $msgpack
done

executing C4H10O_6b12cd9fa9bd56c0692fcd6384bbc05af96c1c1e_e746222796fc2c4c5a1f896fa1cc1cefffe7044c.msgpack
executing C4H10O_8f405ba6ddee36ab6dcc30af38dfc4fea14c5a86_e746222796fc2c4c5a1f896fa1cc1cefffe7044c.msgpack

Now if we run the job again, it still exits.

[15]:

job.run()

generate-conformers: 100%|██████████████████████| 1/1 [00:00<00:00, 6061.13it/s]

An exception has occurred, use %tb to see the full traceback.

SystemExit: Exiting to allow running QM computations; commands are in no_server_working_directory/single_point/run_single_point.sh

/Users/lily/anaconda3/envs/tmpenv3.9/lib/python3.9/site-packages/IPython/core/interactiveshell.py:3377: UserWarning: To exit: use 'exit', 'quit', or Ctrl-D.
  warn("To exit: use 'exit', 'quit', or Ctrl-D.", stacklevel=1)

However, checking the script that PsiRESP generates, only the remaining, uncomputed calculation is listed.

[16]:

%%bash

cat no_server_working_directory/single_point/run_single_point.sh

#!/usr/bin/env bash
psi4 --qcschema C4H10O_1c5e35d36810e33e916c6644f1838761efd45846_e746222796fc2c4c5a1f896fa1cc1cefffe7044c.msgpack

[17]:

%%bash

cd no_server_working_directory/single_point && bash run_single_point.sh

Final run

Now, if we run the job a final time, the ESP files will be read, and the charges will be computed.

[18]:

charges = job.run()

generate-conformers: 100%|██████████████████████| 1/1 [00:00<00:00, 6797.90it/s]
compute-esp: 100%|████████████████████████████████| 3/3 [00:00<00:00, 19.14it/s]

[19]:

charges

[19]:

[array([-0.0629820947554362,  0.1148235481698403,  0.0295351606653119,
         0.1564346120814787, -0.5986354026115893,  0.0071470774475118,
         0.0071470774475118,  0.0071470774475118, -0.0285110420516355,
        -0.0285110420516355, -0.002700823115586 , -0.002700823115586 ,
         0.0032398644142937,  0.0032398644142937,  0.3953269456137151])]

[ ]: