July 2022¶
Friday 1st July¶
I’ve finally managed to get myself vaguely functional, but plonking myself at a bar at 7pm on a Friday night and trying to do a bit of coding. I decided to bit of “play”, and wrote up a little psithon script to compare the effects of certain parameters on the emission of BODIPY.
I’m comparing UFF, Bondi and Alinger radii sets with 1.2x scaling on and off, and am also testing the effect of the GEPOL tesserae area on the excitation energy and PCM polarisation energy.
| Radii Set | 1.2x Scaling | Method | Emission (eV) | Error (eV) | PCM polarisation energy (kj/mol) |
|---|---|---|---|---|---|
| Experiment | 2.407 | ||||
| Bondi | True | IEFPCM | 2.570 | 0.163 | -47.18 |
| Bondi | True | CPCM | 2.563 | 0.156 | -48.36 |
| Bondi | False | IEFPCM | 2.422 | 0.015 | -75.69 |
| Bondi | False | CPCM | 2.414 | 0.007 | -77.24 |
| UFF | True | IEFPCM | 2.660 | 0.253 | -31.34 |
| UFF | True | CPCM | 2.655 | 0.248 | -32.06 |
| UFF | False | IEFPCM | 2.531 | 0.124 | -50.82 |
| UFF | False | CPCM | 2.524 | 0.117 | -51.91 |
| Allinger | True | IEFPCM | 2.577 | 0.170 | -42.93 |
| Allinger | True | CPCM | 2.571 | 0.164 | -43.69 |
| Allinger | False | IEFPCM | 2.427 | 0.020 | -65.65 |
| Allinger | False | CPCM | 2.418 | 0.011 | -67.04 |
| GEPOL Area | |||||
| 0.05 | 2.637 | 0.230 | 27.67 | ||
| 0.1 | 2.411 | 0.004 | 77.53 | ||
| 0.3 | <= Default for above | 2.414 | 0.007 | 77.24 | |
| 0.5 | 2.418 | 0.011 | 76.45 | ||
| 0.7 | 2.418 | 0.011 | 76.43 | ||
| 0.9 | 2.423 | 0.016 | 76.87 | ||
| 1.1 | 2.421 | 0.014 | 77.28 | ||
| 1.3 | 2.421 | 0.014 | 77.28 | ||
| 1.5 | 2.421 | 0.014 | 77.28 |
The results seem pretty straightforward, Bondi radii without scaling, using CPCM seems to be the best for BODIPY, with the smallest GEPOL point area result giving the most accurate results.
In re-doing the area test with the best radii set combo, the GEPOL area is a bit more complicated than before, with the most fine GEPOL area resulting in much worse results.
Tuesday 5th July¶
Decided to give the same process a go with the typical nile red dye as well (emission 657nm in water)
The results are interesting for sure… The optimisation process in orca using gaussian points gave transition energies closer to 2.2 eV. I’m going to o a TD-DFT SP trying to replicate the optimal configuration I found in Psi4.
| Radii Set | 1.2x Scaling | Method | \(k\) | Emission (eV) | Error (eV) |
|---|---|---|---|---|---|
| Experiment | 1.89 | ||||
| Bondi | False | IEFPCM | 0 | 1.84018337 | 0.04981663 |
| Bondi | False | IEFPCM | 0.5 | 1.84018337 | 0.04981663 |
| Bondi | False | CPCM | 0 | 1.8691569 | 0.0208431 |
| Bondi | False | CPCM | 0.5 | 1.87820738 | 0.01179262 |
| Bondi | True | IEFPCM | 0 | 2.01128879 | 0.12128879 |
| Bondi | True | IEFPCM | 0.5 | 2.01128879 | 0.12128879 |
| Bondi | True | CPCM | 0 | 2.01184468 | 0.12184468 |
| Bondi | True | CPCM | 0.5 | 2.01720997 | 0.12720997 |
| UFF | False | IEFPCM | 0 | 1.97683004 | 0.08683004 |
| UFF | False | IEFPCM | 0.5 | 1.97683004 | 0.08683004 |
| UFF | False | CPCM | 0 | 2.01030428 | 0.12030428 |
| UFF | False | CPCM | 0.5 | 2.01581121 | 0.12581121 |
| UFF | True | IEFPCM | 0 | 2.12364278 | 0.23364278 |
| UFF | True | IEFPCM | 0.5 | 2.12364278 | 0.23364278 |
| UFF | True | CPCM | 0 | 2.12175772 | 0.23175772 |
| UFF | True | CPCM | 0.5 | 2.12581901 | 0.23581901 |
| Allinger | False | IEFPCM | 0 | 1.86682167 | 0.02317833 |
| Allinger | False | IEFPCM | 0.5 | 1.86682167 | 0.02317833 |
| Allinger | False | CPCM | 0 | 1.90906567 | 0.01906567 |
| Allinger | False | CPCM | 0.5 | 1.91597263 | 0.02597263 |
| Allinger | True | IEFPCM | 0 | 2.03413434 | 0.14413434 |
| Allinger | True | IEFPCM | 0.5 | 2.03413434 | 0.14413434 |
| Allinger | True | CPCM | 0 | 2.04665886 | 0.15665886 |
| Allinger | True | CPCM | 0.5 | 2.05164165 | 0.16164165 |
| GEPOL Area | |||||
| 0.05 | 1.92563237 | 0.03563237 | |||
| 0.1 | 1.90499419 | 0.01499419 | |||
| 0.3 | <= Default for above | 1.90816755 | 0.01816755 | ||
| 0.5 | 1.90840371 | 0.01840371 | |||
| 0.7 | 1.90163936 | 0.01163936 | |||
| 0.9 | 2.2734716 | 0.3834716 | |||
| 1.1 | 2.2734716 | 0.3834716 | |||
| 1.3 | 2.2734716 | 0.3834716 | |||
| 1.5 | 2.2734716 | 0.3834716 |
In taking the optimised geometry and manually inputting in the Bondi radii for the elements, with and selecting gepol_ses for the cavity, with the COSMO \(\varepsilon\) function.
I tried to recreate the optimal configuration in ORCA, but the results were still not great (\(E_{S^1\to S^0}=2.041\:eV\)) with an error of 0.15 eV. This persists even when making sure that all the permitvities, radii, cavities, etc are the same.
In double checking, the difference between Psi4 and ORCA for gas phase is only ~0.001eV, so the errors have to be occurring in the PCM itself.
TDA seems to be causing errors here. I’ve turned it off and now the results are more in line with each other. I’m testing with the Gaussian point charge scheme in ORCA again to see if the discrepancies return, but the errors I was getting when using the same PCM parameters were within 0.02 eV.
The default vdw_gaussian scheme gives errors as well, of ~0.12eV even with TDA off.
Thursday 7th July¶
I did go into the office today and do a bit of reading and writing, but only for a few hours. I’m definitely starting to get a better feel for PCM and vertical excitations though, which is a good feeling!… though It’s mostly limited to TD-DFT and a bit of CASSCF atm.
Friday 8th July¶
In thinking more closely about the radii problem above, it seems really tricky to break down the problem because the PCM method on its own only accounts for the electrostatic contribution, and misses out on repulsion-dispersion and cavitation. Any optimisation of the electrostatic term to match experiment will be overestimating electrostatics to compensate for the other terms.
Should we just assume that the n-el terms cancel out in the vertical excitation? the accuracy of SMD seems to suggest that we can’t.
Monday 11th July¶
I noticed as I was explaining to a friend what a Hessian is that some of my radicals are not at their minimum structures, so I’ve gone through, grep’d out those files and have resubmitted an opt that reads in that hessian.
Just found FPBase which is a fantastic absorbance and emission spectra viewer for some prototypical chromophores and fluorophores.
Tuesday 12th July¶
I’ve spent my day trying to build the dataset. I’ve found a few florurophores that are convenient and well documented and I’ve been putting them into a csv file.
I wrote a script to take the smiles strings and convert them to 3D and run a conformer search.
Wednesday 13th July¶
I changed up my test cases a little bit and found reference data for all the species I’m going to check
I did a CREST search on all of these structures and have queued up a bunch of SMD ground state optimisations with vertical excitations. When that’s done, I’ll get the excited state geometries and excitations as well, then the fun can begin!
I should probably get back to reading…
Thursday 14th July¶
I spent my day plonked at ACMI and just read reviews. I’ve been working my way through Dielectric continuum methods for quantum chemistry. It’s a pretty good review as it seems to break things down and step through the maths in a very palatable manner.
I had a little bit of a play with trying to figure out how to set up my data and code structure for my PhD, which I still need to figure out how to do well. At the moment, everything is just stored in a CSV file, that I can import into pandas, but that’s going to get messy quite quickly.
Friday 15th July¶
Today I’ve not really done any work. I ran maintenance on the jobs on MonARCH to make sure everything is still running, but apart from that, I’ve been running errands all day.
Tonight I’m going to a [public
with a couple of friends, at Melbourne Uni on quantum computing, which should be fun 
Monday 18th July¶
I’ve been working on transferring my work to an ipynb workflow that involves being able to pull and push data from the cluster remotely and process everything in a local pandas dataframe. I’ve also been functionalising a lot of this code to make it more extensible. Today is clearly a code-heavy day, which I guess means I should focus on reading tomorrow.
I also had some training today for CHM1022.
Wednesday 20th July¶
Today and yesterday I mostly focused on code and job management, and started building a module to be able to interface cleanly with MonARCH via SSH without spawning a huge number of sessions (and getting blocked).
I’ve fired off my benchmark jobs again, this time including bot dz and tz basis sets
Thursday 21st July¶
Today I’m back at ACMI. I’ve picked hi01 as the node that I’m going to license QChem to for testing purposes, so I’ve just sent that off. I found a few errors in my code that needed fixing, but I should really just focus on reading now…
Friday 22nd July¶
Today I’ve been playing with Q-Chem now that I’ve got the install on hi01 running and trying to get some jobs running. It’s not as black box as ORCA, so I’m having to manually specify grids and RI parameters and things.
I’ve figured out a few things, like:
rpa = 0is TDA,rpa = 1is TDA as a guess for RPA, andrpa = 2is pure RPA- It will automatically calculate triplets (which will change what root 1 is) unless you specify
CIS_TRIPLETS false - Density fitting seems to be automatically enabled when you set an
AUX_BASIS - Memory setting is optional but is set with
MEM_TOTALand is specified in MB for the whole job GUI 2specifies to output an IQMol readable.fchkfile.- Integration gids can either be set explicitly with
XC_GRID 000099000590, or can be set to a preset with the grids in the table below
| Grid | Radial/spherical | Default for functionals |
|---|---|---|
| SG-0 | (23, 170) | None |
| SG-1 | (50, 194) | LDA, most GGAs and hybrids |
| SG-2 | (75, 302) | Meta-GGAs; B95- and B97-based functionals |
| SG-3 | (99, 590) | Minnesota functionals |
I also updated my 2slm script to be able to handle qchem input files.
Monday-Wednesday 25-27th July¶
I’ve spent the past few days reading, pushing through jobs still and trying to get some baseline data. The REALLY tricky part of this dataset work seems likely to be that the geometry for each species is going to be solvation dependent, so using the one optimised geometry is likely not going to work.
Even if I were to optimise the molecules at EOM-CCSD or ADC(2) or something fancy and expensive, the solvation and excited state handling of solvation would make it impossible to have a “gold standard” option.
I had a look at Q-Chem EFP, and it looks like you can only do:
- EFP only optimisation (no QM component at al)
- EFP and QM SP
This mean that I couldn’t use the iterative xtb/efp method that I developed for use with Psi4, since I can’t get a QM gradient with EFP fragments. Maybe I need to fire up GAMESS for that? :/
Friday 29th July¶
I had a good meeting with Katya (see notes below) and in the afternoon I helped a third year PhD (Sarmi) that wanted to be able to run some calculations to validate their work.
Saturday 30th July¶
I’m plonked in ACMI and have diced to play with some QMMM in Q-chem, which means atomtyping my molecules. I’m using ligpargen to get the OPLS Parameters
Important (play) coding:¶
- Write a server-side script for reading files? (should be faster than catting over ssh)
- Potentially explore sftp file reads instead of ssh catting
- QM/MM in Q-Chem
- QM solvation shell + EFP for long range
- Having a few explicit water molecules around charged functionality + EFP/PCM
From meeting with Katya:¶
- Perhaps we need to do an assessment study comparing basis set types with associated errors. Pople basis sets do seem to be more commonly used
- Perhaps we could use a smaller integration grid as well? (75/302)
- Probably better for CAM-B3LYP than ωB97X-D3
- Perhaps we could use a smaller integration grid as well? (75/302)
- ML might be the way to go to get non-linear regressions
- RT-TDDFT might be a further consideration if that’s how I want this PhD to go, but probably isn’t worth it for now.
To ask Toby:¶
- Is my dataset complete enough?
- Is there more space of fluorophores we need to consider?
- Would it be possible to get our own spectrochemical data?
To do:¶
Apart from all the papers that need reading…
- Finish the Benedetta Mennucci review
- Finish the non-Benedetta Mennucci review
- Read up on domain decomposition solvation (ddPCM/ddCOSMO)
- Read up on COSMO-RS
-
Should we get a license? Perhaps not, since this would be effectively paywall restricted science…- Yeah, not worth it. Katya agrees.
-
- Read program’s manual’s solvation sections and pull out key discussed papers
- ORCA
- Q-Chem
- PCMSolver (read the paper)