Visualize canonical orbitals from CPMD¶
For different reasons, one might be interested in looking at the Kohn-Sham canonical orbitals. In my case it was to identify the main orbitals involved TDDFT electronic transitions.
Step 1: Get the CPMD files for the orbitals to plot¶
The CPMD keyword used to generate the required file is RHOOUT BANDS in the &CPMD section.
I have written a script that reads a CPMD TDDFT output file
and extract information regarding the orbital involved in the
TDDFT electronic transitions.
The script is called read_orbs_info
and is part of the comp_chem_py package.
The script prints some of the required lines to plot orbitals with CPMD.
For example:
RESTART WAVEFUNCTION COORDINATES LINRES LATEST
KOHN-SHAM ENERGIES
20
RHOOUT BANDS
4
-138 -139 -140 -141
Those lines have to be pasted in the &CPMD
section of a CPMD input file.
The restart information and number of Kohn-Sham energies is set to match
the TDDFT output file that is read by the script.
Finally, the lines after RHOOUT BANDS
are set such that a grid file will
be generated for all the orbitals involved in the TDDFT electronic transitions.
Of course all the information can be modified to suit your purpose.
Step 2: Generate cube files¶
After the first step we should have a bunch of WAVEFUNCTION.index files available. 8 in our case. Those files have to be converted to the Gaussian cube format. This can be done using the cpmd2cube.x utility. This script is provided with the CPMD package and is located at:
/share/lcbcpc35/software/cpmd2cube/cpmd2cube.x
On the LCBC cluster.
To use it just run:
cpmd2cube.x WAVEFUNCTION.*
Or use the -h
option for help.
This operation should produce a cube file and a pdb file for each input, e.g.
WAVEFUNCTION.8
leads to WAVEFUNCTION.8.cube
and WAVEFUNCTION.8.pdb
.
The cube files can then be visualized using e.g. the VMD package:
module load vmd
vmd WAVEFUNCTION.8.pdb
Then in the VMD terminal window:
vmd > mol addfile WAVEFUNCTION.8.cube
vmd > mol addrep 0
vmd > mol modmaterial 1 0 Opaque
vmd > mol modstyle 1 0 Isosurface 0.01 0 0 0 1 1
vmd > mol modcolor 1 0 ColorID 0
vmd > mol addrep 0
vmd > mol modstyle 2 0 Isosurface -0.01 0 0 0 1 1
Or using the vmd_plot_cube
script from the comp_chem_py package:
vmd_plot_cube GEOMETRY.xyz -c WAVEFUNCTION.*.cube