Plotting wave functions

Creating a wave function file

The following script will do a calculation for a CO molecule and save the wave functions in a file (CO.gpw).

from ase import Atoms
from gpaw import GPAW

d = 1.1   # bondlength of hydrogen molecule
a = 5.0   # sidelength of unit cell
c = a / 2
atoms = Atoms('CO',
              positions=[(c - d / 2, c, c),
                         (c + d / 2, c, c)],
              cell=(a, a, a))

calc = GPAW(mode='fd', nbands=5, h=0.2, txt=None)
atoms.calc = calc

# Start a calculation:
energy = atoms.get_potential_energy()

# Save wave functions:
calc.write('CO.gpw', mode='all')

Creating wave function cube files

You can get separate cube files (the format used by Gaussian) for each wavefunction with the script:

from ase.io import write
from ase.units import Bohr
from gpaw import restart

basename = 'CO'

# load binary file and get calculator
atoms, calc = restart(basename + '.gpw')

# loop over all wfs and write their cube files
nbands = calc.get_number_of_bands()
for band in range(nbands):
    wf = calc.get_pseudo_wave_function(band=band)
    fname = f'{basename}_{band}.cube'
    print('writing wf', band, 'to file', fname)
    write(fname, atoms, data=wf * Bohr**1.5)

The script produced the files CO_0.cube .. CO_5.cube, which might be viewed using for example jmol or VMD.

Plotting wave functions with jmol

See http://jmol.sourceforge.net/docs/surface/

You can visualize an isosurface from a cube file by creating the \(myscript.spt\) script:

load CO_1.cube

#isosurface name cutoff value cubefile
isosurface pos cutoff 0.05 CO_1.cube
color isosurface translucent red
isosurface neg cutoff -0.05 CO_1.cube
color isosurface translucent blue

and executing it with jmol:

jmol -s myscript.spt

You can also save the image directly on the command line:

jmol -ions myscript.spt -w JPEG:myfile.jpg

Plotting wave functions with VMD

To view the wavefunctions, start VMD with the command line:

vmd CO_*.cube

You will get two windows, one with a very primitive representation of the molecule, and one called ‘VMD Main’. In the VMD Main window, select Graphics/Representation. In the top panel there will be a single “representation” of the molecule called “Line”. Select it, and change the type to “CPK”. This will show the molecule with spheres for atoms and rods for bonds. You probably want to reduce the size of the spheres and increase the resolution of the spheres. Now you see the atoms!

To see the wavefunctions, click on the “Create Rep” button, and select the new representation. Select type “Isosurface”. Near the bottom, find the menu labeled “Draw” and select the value “Solid Surface”. Now you can see an iso-surface of the wavefunction, you select the value in the field “Isovalue”. The default 0 is rarely useful. You select the different Kohn-Sham states (stored in the different cube files) with the pull-down menu labeled “Vol”.

IMPORTANT: This works best for molecules. In solids, the wavefunctions will be complex, VMD does not handle this well.