Plotting charge densities

Command Summary

The following commands are sufficient to run the tutorial.

Preliminaries

This tutorial assumes you have cloned and built the lmf repository (located here). For the purpose of demonstration, ~/lmf will refer to the location of the cloned repository (source directory). In practice, this directory can be named differently. It is assumed you have moved the binaries to a directory in your path.

This test case demonstrations generation of charge density in a plane:

$lmf/b/./fp/test/test.fp --quiet --poszer --MPIK bzt  This tutorial goes through that demonstration, with some additional explanation. It assumes you will run in MPI mode, and uses prefix mpirun -n 9 for MPI jobs. To run in scalar mode, you can skip the prefix. At present, there is no capability to interpolate the smoothed density to an arbitrary plane, so you are restricted to choosing a plane that passes through uniform mesh points. (The interstitial density is stored on a uniform grid). To make this translation a little easier, you can run lmf with switch --wden~atpos. This doesn’t generate a density but prints out a table indicating where the nuclei are as multiples of the mesh points. Tutorial Start in a fresh working directory of your choice (better not use use the directory where you compiled the code). Copy the input files needed to run this calculation from the distribution: $ cp ~/lmf/fp/test/{ctrl.bzt,site.bzt} .


Make the free atom density and the self-consistent density

 $lmfa bzt$ mpirun -n 9 lmf bzt -vnit=30


Get information about the site positions as multiples of the uniform mesh

$lmf bzt --rs=0 -vnit=0 -wden:atpos  You should see a table like this:  Site positions as (possibly fractional) indices of 36x36x44 smooth mesh Site p1 p2 p3 Site p1 p2 p3 Site p1 p2 p3 1 0.0 0.0 0.0 | 6 0.0 0.0 22.0 | 11 29.8 23.7 14.3 2 12.0 24.0 14.9 | 7 0.0 18.0 0.0 | 12 29.8 6.2 14.3 3 24.0 12.0 29.1 | 8 18.0 18.0 0.0 | 13 12.3 6.2 14.3 4 12.0 24.0 36.3 | 9 18.0 0.0 0.0 | 14 23.7 29.8 29.7 5 24.0 12.0 7.7 | 10 6.2 12.3 29.7 | 15 6.2 29.8 29.7  Note the lattice vectors are hexagonal, as can be seen from out put of, e.g., lmchk  Plat Qlat 0.000000 -1.000000 0.000000 0.577350 -1.000000 0.000000 0.866025 0.500000 0.000000 1.154701 0.000000 0.000000 0.000000 0.000000 1.227430 0.000000 0.000000 0.814710 alat = 10.926401 Cell vol = 1386.624179  The following writes writes out the density generated by Mattheis construction (superposition of free atomic densities) in the basal plane and passing through the origin, to file atrho.bzt. $ lmf bzt --rs=0 -vnit=0 --wden:core=2:l1=1,0,0,37:l2=0,1,0,37:fn=atrho


The syntax of --wden is described here; in this context they mean the following:

• core=2    tells lmf to exclude core densities.
• ~l1=1,0,0,37  specifies the first axis for the 2D mesh (1,0,0) corresponds to (0,-1,0) in Cartesian coordinates (see Plat above).
It selects 37 points because the mesh along P1 consists of 36 points. Thus points 1 and 37 are equivalent.
• ~l1=0,1,0,37  specifies the second axis for the 2D mesh, which is $(\sqrt{3}/2,1/2,0)$. 37 points are specified because P2 also has 36 points.
• fn=atrho   specifies the file name for file I/O.

When writing file atrho.bzt, lmf prints the following:

 ioden : local densities + envelope density, Qloc=19.189235  Q=111.999997
Writing smooth density to file atrho : origin at (0,0,0).
Origin, in cartesian coordinates  0.000000 0.000000 0.000000
v1: (36+1 pts) = (1,0,0)(p1,p2,p3) = (0.000000,-1.000000,0.000000) l=1.000000
v2: (36+1 pts) = (0,1,0)(p1,p2,p3) = (0.866025,0.500000,0.000000) l=1.000000
Angle between vectors : 2.094395 radians = 120 deg

Unit normal vector (0.000000,0.000000,1.000000).  Origin lines in plane 0
NN planes connected by (0,0,1)(p1,p2,p3) = (0.000000,0.000000,0.027896)
Planes separated by 0.0278961, 44 planes to shortest period, 0 planes in normal direction


v1 and v2 are the axes for the plane: they are displayed as multiples of Plat and also in Cartesian coordinates. Note the angle between v1 and v2 is 120o.

Next we write the self-consistent density in the same plane, this time to the default file nane smrho.bzt.

$lmf bzt --wden:core=2:l1=1,0,0,37:l2=0,1,0,37  Finally let’s draw a contour plot of the charge density. Rather than plot the density itself, we can plot the shift in the density relative to the Mattheis construction: this shows where the bonding goes. $ mcx smrho.bzt atrho.bzt -- > rho_diff.bzt


To draw a figure, copy the following box to file plot.rho

fplot

# ... for (sqrt(3)/2,1/2,0), (0,-1,0)
-frme:theta=2*pi/3 .4,1.3,0,.9

-lt 1,bold=3,col=.00,.1,1.,fill=0 -con  .001 -qr rho_diff.bzt
-lt 1,bold=3,col=.25,.1,.7,fill=0 -con  .002 -qr rho_diff.bzt
-lt 1,bold=3,col=.50,.1,.4,fill=0 -con  .005 -qr rho_diff.bzt
-lt 1,bold=3,col=.75,.1,.1,fill=0 -con  .010 -qr rho_diff.bzt

-lt 3,bold=5,col=.1,.00,1.0,fill=0 -con -.001 -qr rho_diff.bzt
-lt 3,bold=5,col=.1,.25,.75,fill=0 -con -.002 -qr rho_diff.bzt
-lt 3,bold=5,col=.1,.50,.50,fill=0 -con -.005 -qr rho_diff.bzt
-lt 3,bold=5,col=.1,.75,.25,fill=0 -con -.010 -qr rho_diff.bzt
-lt 3,bold=5,col=.1,1.0,.00,fill=0 -con -.020 -qr rho_diff.bzt


Then run

$fplot -f plot.rho$ open fplot.ps


(Use your favorite postscript viewer for open).

Solid contours correspond to positive densities (actually density differences), varying from blue for most intense, to red for least intense. while dashed contours are for negative densities, varying from blue to green. It is apparent that charge is piling up in a rim around the O atom, and being taken away from a rim around the Ba, consistent with charge transfer from cation to anion in a polar compound. Click below to view the figure fplot should generate.

In units of the lattice vectors, the top right corner is (1,1,0). Atom 8 (an O atom) is located at (1/2,1/2,0) in these units, which lies at the centre of the figure. A Ba sits at the bottom left corner. You can readily see this by invoking

$lmchk bzt --pr50  You should see  site spec pos (Cartesian coordinates) pos (multiples of plat) 1 Ba2 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 ... 8 O1 -0.433013 0.250000 0.000000 -0.500000 -0.500000 0.000000  Plotting the effective potential You can, instead of plotting the density, plot the (pseudo) effective potential. The true one has singularities at each nucleus; to smooth it out, an effective one is written. lmf bzt --rs=0 --wden:core=2:l1=1,0,0,37:l2=0,1,0,37:fn=atrho:veff lmf bzt --wden:core=2:l1=1,0,0,37:l2=0,1,0,37:veff mcx smrho.bzt atrho.bzt -- > rho_diff.bzt  To draw a picture, cut and paste the box below to file plot.pot fplot # ... for (sqrt(3)/2,1/2,0), (0,-1,0) -frme:theta=2*pi/3 .4,1.3,0,.9 -lt 1,bold=3,col=.00,.1,1.,fill=0 -con .01 -qr rho_diff.bzt -lt 1,bold=3,col=.20,.1,1.,fill=0 -con .02 -qr rho_diff.bzt -lt 1,bold=3,col=.40,.1,1.,fill=0 -con .05 -qr rho_diff.bzt -lt 1,bold=3,col=.60,.1,1.,fill=0 -con .1 -qr rho_diff.bzt -lt 1,bold=3,col=.80,.1,.7,fill=0 -con .2 -qr rho_diff.bzt -lt 1,bold=3,col=.99,.1,.4,fill=0 -con .5 -qr rho_diff.bzt -lt 3,bold=5,col=.1,.00,1.0,fill=0 -con -.01 -qr rho_diff.bzt -lt 3,bold=5,col=.1,.25,1.0,fill=0 -con -.02 -qr rho_diff.bzt -lt 3,bold=5,col=.1,.50,1.0,fill=0 -con -.05 -qr rho_diff.bzt  Then run $ fplot -f plot.pot
\$ open fplot.ps