c2x and Molecules

One of the better programs for plotting isosurfaces appears to be XCrysDen. It is not perfect, but it tends to do what I need. There are plenty of alternatives.

Let us start with ethene in an 8x8x7 Angstrom box. Sufficient input files for CASTEP are:

8.000000 8.000000 7.000000
90.000000 90.000000 90.000000
%endblock LATTICE_ABC

  C 4.698560 4.000000 3.500000
  C 3.301440 4.000000 3.500000
  H 5.270880 4.923920 3.500000
  H 2.729120 3.076080 3.500000
  H 2.729120 4.923920 3.500000
  H 5.270880 3.076080 3.500000

kpoints_mp_grid 1 1 1
kpoints_mp_offset 0.25 0.25 0.25
cut_off_energy  = 40 ry
iprint          = 1
opt_strategy : speed

Castep should generate its own pseudopotentials, and run that calculation in about a minute. Note that the k-point chosen will minimise interactions between the neighbouring periodic images.

Charge Density

$ c2x -v -c ethene.check ethene.rho.xsf
$ xcrysden --xsf ethene.rho.xsf &
Ethene ball and stick

Try clicking the "Zoom +" button a few times, then select Tools, Data Grid. The input file contains just a single data grid, which will be selected, so click "OK", type in an isovalue of 0.5 in the next dialogue box, and then click "Submit". Drag with the mouse to rotate, and one bone, or ethene charge density, results.

Ethene charge density

Should you feel tempted to save or print the image, you may wish to change the background to white. Click on the coloured box at the top left beside the "File" menu, and select white.

Unfortunately XCrysDen limits cut planes to being parallel to the cell axes, but they can still look pretty. In this case the bonds have been removed with Modify, Ball/Stick Ratio, and setting the ratio to zero.

Ethene charge cut plane


Of course ethene has a famously pretty pi bond. A little basic chemistry suggests that this will be the highest energy band, and there are clearly six bands in total, so it can be extracted with:

$ c2x -v -b=6 ethene.check ethene.pi.xsf
$ xcrysden --xsf ethene.pi.xsf &

(Check2xsf will complain about a few thousand imaginary components. As long as this is very much less than the size of the FFT grid, it probably isn't important. Using a larger box, and tigher electronic minimisation convergence criteria, would reduce this number.)

Ethene pi bond

For the sigma bonds, it is possible to create a single xsf file containing all the sigma bands, using -b=1-5. Unfortunately, XCrysDen insists on rewriting its input file on start-up, and rewriting this file can take several seconds. It does then allow one to choose any linear combination of the bands to plot, but one cannot change one's choice without restarting XCrysDen.

A naive chemist might assume that there was one C-C sigma bond, and four degenerate C-H sigma bonds. This is not so: these bonds will interact, and the sigma-like bands found by CASTEP are as shown below (after stitching together with the Gimp):

Ethene sigma MOs

It is possible to extract eigenvalues via check2xsf. It needs a relatively high verbosity, and it needs to be prompted into looking at the wavefunctions.

$ c2x -vv -b=1 ethene.check /dev/null
is sufficient.

If one wishes to view charge densities arising from a band, rather than psi, then specifying B rather than b should achieve this. One can also cause check2xsf to sum bands, rather than relying on the visualisation software to do so.

-B=2,3,5 -A
will sum the densities from bands 2,3 and 5. It is not possible to specify prefactors.