back to summary
-
You can access the necessary files for this tutorial by downloading this
archive. Once downloaded, place it in a convenient location and extract its contents:
[user]$ tar -xvjf 4_intramol_piStack.tbz2
Unpacking the archive will create a directory called 4_intramol_piStack/, wherein this tutorial will be conducted. - The results of the calculations can be found in this archive_results.
IV. Quantifying intramolecular π-stacking: QM mode (wfx/wfn)
Go to directory 4_intramol_piStack/We are now going to study the intramolecular interaction between two aromatic rings in a single molecule. Prior to the calculations, you may have a look at the molecule below or load it in VMD from the available .xyz file.
In the QM mode, the simultaneous combination of FRAG1 and FRAG2 keywords in param.igm allows to probe interactions between two subsets of a given molecule. In order to do so, both FRAG1 and FRAG2 atom selections must be supplied in the param.igm input file.
Edit and save the following param.igm input file:
1
mol.wfx
FRAG1 11-21;35-42;
FRAG2 1-6;26-30;
FRAG1 11-21;35-42; (oxindole)
FRAG2 1-6;26-30; (phenyl)
- 1 --> one single file describes the system
- mol.wfx --> the name of this file --> here, the system is described at the QM level of theory
(when the .wfx or .wfn file format is detected)
Provided that you have installed IGMPlot (see documentation for installation instructions): run the program IGMPlot with the command in a linux terminal: -> IGMPLOT param.igm > igm.log &
- IGMPLOT: the executable name of the IGMPlot program
- param.igm: the input file
- > igm.log: standard output will be redirected to file igm.log
- &: command is run in background
This calculation will take time (around 6 minutes on a 4-CPU processor). Have a look at the “runinfo” information file giving you an estimate of the remaining time before the end of calculation:
Have a look at the resulting output file list in the current directory:
- AtContribInter.vmd: vmd script to display the atoms coloured according to their contribution to the intermolecular interaction
- mol-AtContribInter.dat: atomic contributions to the intermolecular interaction
- igm.log: text outputfile of the IGM execution
- igm.vmd: VMD script to automatically display δg 3D isosurfaces
- mol-dens.cube: cube file describing the electron density (ED) at each node of the grid
- mol-dgInter.cube: cube file describing the δg_inter value at each node of the grid
- mol-dgIntra.cube: cube file describing the δg_intra value at each node of the grid
- mol-igm.dat: file with data in columns to be used with GNUPLOT to generate 2D-fingerprints
- mol-nci.dat: file with data in columns to display the RDG=f(ED) 2D-plot
- mol-RDG.cube: cube file describing the Reduced Density Gradient value at each node of the grid
- nci.vmd: VMD script to automatically display RDG 3D isosurfaces
3D-isosurfaces
We are going now to use the igm.vmd file generated by IGMPlot to identify the regions of space where interactions take place between the two intramolecular fragments. This vmd script must be located in the folder containing the cube files while running VMD. Run VMD and load the igm.vmd session file (menu File/load Visualization state). Two windows appear: the VMD main window and the molecular viewer.
In the VMD main panel appear two “Molecules” lines (res-dgInter.cube and res-dgIntra.cube). Default of the script igm.vmd is to display the interactions between fragments: resdgInter.cube Molecule (button D is activated). Here, default is then to display the interactions between two fragments.
The molecular system appears in the molecular viewer in “ball & sticks” representation. Press “S” to scale the picture, “T” to translate or “R” to rotate.
As you can see (on the left on the Figure), the flat region of interaction between the two aromatic rings could be extracted from the wave function content.
Since δg captures ED clash between different sources (atoms, fragments, molecules) such δg isosurfaces give rise to a coherent picture with stabilizing interactions accumulated in the center of several envelopes.
Note: the default δg isovalue = 40% of the maximum peak height found by IGMPlot in the 2D-finger-print. Changing the iso-value will change the size of the isosurfaces. To do that, go to the menu Graphics/Representations, select the “Isosurface” representation (by clicking on it) and change the “Isovalue” number in the panel, to gradually grow or diminish the volume enclosed in the envelops. The larger the isovalue, the larger the ED interaction.
♥ Always check the intra-representation to make sure that your FRAG1/FRAG2 definition meets your needs. To do so, double-click in the VMD main menu on the "D" letter facing the res-dgIntra.cube system, and simultaneously disable the res-dgInter.cube system (by double-clicking in the VMD main menu on the "D" letter facing the res-dgInter.cube system). This yields the δg_intra 3D-isosurfaces represented at the center of the Figure above.
Finally, have a look at the atomic contributions to the interactions between the two aromatic rings (on a blue to red color scale): exit VMD, re-run VMD and load the script AtContribInter.vmd (representation on the right).
Note that, very rarely (we observed it only once for a radical TS species), the use of diffuse functions in the basis set might generate spurious δg_inter isosurfaces. In such case, simply change the basis set (setting no diffuse function), the IGM scoring being hardly influenced by the basis-set in general cases. An alternative is to use in this case the hybrid QM/Promolecular Hirshfeld approach (see documentation, HIRSH keyword).
Analysis of igm.log
The corresponding output emphasizes how many orbitals have participated to probe this internal interaction: here, only those primitives belonging to the two fragments have been used:
Next, an "INTER" section is reported to describe the interaction detected between the two fragments. Note immediately that IGMPlot has detected the MAXIMUM of the δg peak(s) = "2.3447e-02 a.u.", classifying this interaction in the “non-cov” domain on the indicative IGM peak scale:
Noteworthy, the magnitude of this π-stacking interaction lies in between that of the typical hydrogen-bond in the water dimer (δg_peak = 0.06 a.u.) and that of the π-stacking present in the benzene dimer (δg_peak = 0.01 a.u.).
Next, a number of integrated scores are reported. Herafter are described the meanings of the scores given in the igm.log output.
The score [1] corresponds to the full integration of δg. Thus, every point in the grid is taken into account, and all the δg values are summed and weighted by the elementary volume dv of each grid cell. The [1a] score is the same as [1] except that only those points of the grid for which λ2 > 0 are considered. The [1b] score is the same as [1] except that only those points of the grid for which λ2 < 0 are considered. Thus, [1] = [1a] + [1b]. The separated scores [1a] and [1b] have been implemented upon user request since the sign of λ2 is used in the original NCI paper by Yang to characterise the kind of interaction in the interatomic regions (bonding situation when λ2 < 0 and non-bonding when λ2 > 0).
Next, the score [2] (similar to [1b]) corresponds to the integration of grid points for which λ2 < 0 and for which the signed electron density rho : -0.09 < ρ < 0.09 and for which qg > 1.3. This score aims at extracting and quantifying only weak interactions in a given 2D-plot signature. Here, IGMPlot detects some weak interactions (π-stacking here) between the two fragments. The latter condition allows to get a score focusing only on the peaks of the 2D-plot fingerprint, removing points of low ED contragradience (the flat part of the signal), making the score [2] more sensitive to the interaction situation. The score [3] is similar to [2] but extracts and quantifies strong (covalent) interactions in a given 2D-plot signature. Note that here, absolutely no strong interactions are detected between the fragment (no covalent bonding between the fragments). Of course, the integration of the IGM-δg descriptor cannot be expected to lead as accurate results as QM calculations or IQA, but it offers a fast and robust way to get an order of magnitude of interactions comparing different systems or during their evolution along a reaction. See for instance an example of IGM application in the paper "Mechanistic insights into Smiles rearrangement. Focus on π-π stacking interactions along the radical cascade" (http://doi.org/10.1039/D0OB01511C).