The ability to plot the transmission coefficients in various ways is an important addition to the basic functionality in VNL, so Nori has started a very relevant topic here!
I would like to offer a second version, with focus on a few key points:
The k-point grid. I hereby provide a Monkhorst-Pack class that can be used very generally. The main benefit is that the k-point generation is indepedent of the geometrical structure, so there is no need for complicated reading of the NetCDF file, and moreover there is a well-defined method to integrate over the k-points to obtain the total transmission.
I will not document the Monkhorst-Pack class in this post; see the file, it should be quite self-explanatory.
The plotting. I very much like Nori's plots, 3D plots are always nice. In case you are ok with the contour plot version only, it is however possible to do the plot with ATK itself, thus reducing the need for exporting the data to an external file, and use an external plotting program.
Now, because the Monkhorst-Pack class is so general, the final user script becomes very small. Here is the part which does the calculations:
from ATK.TwoProbe import *
from MonkhorstPack import MonkhorstPackGrid
nc_filename = 'femgofe_antipara.nc'
Npoints = 35
energy = 0.*eV
scf = restoreSelfConsistentCalculation(nc_filename)
kgrid = MonkhorstPackGrid((Npoints,Npoints),False)
K = kgrid.getFractionalKPoints()
trans_coeff = calculateTransmissionCoefficients(scf,energy,(K,Spin.Up))
Modify the code as needed for number of points, the energy, the NetCDF file name, and the spin component (the other one is Spin.Down). If you have an un-polarized calculation, replace the tuple (K,Spin.Up) with just plain K.
Next, we prepare the data for plotting:
Kxy = K[0:Npoints,1]
T = abs(trans_coeff.reshape(Npoints,Npoints))
Finally, the plotting part:
import pylab as P
fig = P.figure()
P.xlabel('kx')
P.ylabel('ky')
P.contourf(Kxy,Kxy,T,40,cmap=P.cm.Spectral)
P.colorbar()
P.axis([-.5,.5,-.5,.5])
P.yticks(numpy.arange(-0.5,0.51,0.1))
P.xticks(numpy.arange(-0.5,0.51,0.1))
P.title('FeMgO, anti-parallel config, majority spin')
P.savefig('femgo_antipara_trans_up.png',dpi=100)
Play around with the plot parameters as you want (title, file name, etc)! The resulting plot of the example above is attached.
(A side-comment: A subtle but hard-to-catch error will occur if you name the Monkhorst-Pack grid variable "grid", which seems like a natural thing to do, and then, instead of import pylab you think, "hey, why not do from pylab import * to make it simpler". The script will then fail, because the pylab import statement overwrites the name "grid" (it's a pylab function). This demonstrates the danger of from X import * statements; we generally recommend the use of import X instead (except I guess for from ATK.KohnSham import * :) !)
The above script is simpler than Nori's (because I have hidden all the dirty stuff in the Monkhorst-Pack class), and therefore a bit easier to use. I think it is also a benefit to use a well-defined Monkhorst-Pack grid, for the reasons mentioned above. Note also that you don't really need the unit cell data to construct the k-point grid for the transmission, since we always use fractional coordinates for that.
The disadvantage is that I did not consider time-reversal symmetry, as Nori did. So, mine is slower by a factor of 2, which probably is not be very critical in many cases, but it might be in some. I will, therefore, show in a separate post how we can take advantage of symmetries in the transmission spectrum, including also mirror planes etc, in a very powerful and time-consuming way.
So, stay tuned - there is more to come!
There is no standardized format for such files, so it depends a bit on which software you want to use. The data you need is kx and ky, which are contained in the variable Kxy (same grid in x and y), and the actual matrix T(kx,ky) which is the variable T.
So, if your plotting software needs a format of the form
x1 y1 z
x1 y2 z
you would loop over Kxy twice and print the matrix:
for x in range(len(Kxy)):
for y in range(len(Kxy)):
print Kxy[x],Kxy[y],T[x,y]
print
The extra print statement gives a new line after each sequence of y values.
If this works or not depends , as mentioned, on the format expected. I guess you just have to play around with it, and try different things (obviously with a small number of k-points at first, not to waste time).
Other software might just need the matrix itself and Kxy (I'm thinking of Matlab for instance). In that case print the matrix T and Kxy as they are, probably to separate files.