Check2xsf and Molecules
Some notes on visualising CASTEP results from aperiodic systems.
One of the best programs for plotting isosurfaces appears to be XCrysDen. I don't think it is perfect, but it tends to do what I need.
Let us start with ethene in an 8x8x7 Angstrom box. Sufficient input files for CASTEP are:ethene.cell
%block LATTICE_ABC ang 8.000000 8.000000 7.000000 90.000000 90.000000 90.000000 %endblock LATTICE_ABC %block POSITIONS_ABS 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 %endblock POSITIONS_ABS kpoints_mp_grid 1 1 1 kpoints_mp_offset 0.25 0.25 0.25ethene.param
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.
$ check2xsf -v -c ethene.check ethene.rho.xsf $ xcrysden --xsf ethene.rho.xsf &
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.
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.
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:
$ check2xsf -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.)
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):
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.
$ check2xsf -vv -b=1 ethene.check /dev/nullis 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 -Awill sum the densities from bands 2,3 and 5. It is not possible to specify prefactors.