First take a look at the analysis objects calculated in the previous section. Locate the file au_dtb_au.nc in the VNL file browser, select the "Transmission Spectrum" object in the Result Browser (upper right panel), and finally click the "Show Plot" button in the lower panel, the Result Viewer.

There are three broad peaks, at around -3, -1.5 and 2.5 eV. Also note the narrow peak at 2.394 eV which rises above . In the following you will investigate the origin of these peaks.

Next select the Device Density of States object in the NC file, and click the Show button to open the DOS analyzer window. Select all carbon and hydrogen atoms using Ctrl left-click. If you zoom into the figure by drawing a rectangle with the left mouse button, you should now see the following plot.

Note the similarity in the peak structure of the PDDOS and the transmission spectrum, indicating that there is a clear correspondence between the energy levels on the phenylene ring and the transmission spectrum.

To understand the energy levels of the phenylene ring in the presence of the surrounding electrodes, you will now calculated the MPSH states. The MPSH states are obtained by diagonalizing the molecular part of the full self-consistent Hamiltonian [1].

Open the Script Generator and

Open the Analysis from File script block and select au_dtb_au.nc from the file browser.

Open the MolecularEnergySpectrum block and select projection on carbon and hydrogen.

Send the script to the Job Manager and execute it (the execution only takes a few seconds). Inspect the MPSH energies in the log file.

Note that state 13, 14, 15 and 16 have energies in the range relevant for the transmission spectrum. In the following you will investigate the symmetry of these states.

To calculate the eigenstates of the MPSH levels, reopen the last Script Generator (it is available from the windows menu at the top of each VNL window, if you didn't close it).

Add four Analysis/EigenState blocks.

Open each EigenState script block and select projection on carbon and hydrogen. Set the Quantum number to 13, 14, 15, 16 in each respective block.

You may delete the MolecularEnergySpectrum block, it is no longer needed (select it and press the Delete key).

Send the script to the Job Manager and execute the script.

Select dtb_analysis.nc in the file browser, and choose the first eigenstate in the result file, and click the button Isosurface>Show. In the Viewer which opens up, open the Plot>Properties menu and select the Isosurface in the Plot tree.

To visualize the eigenstate together with geometry, drag and drop the file au_dtb_au.py onto the Viewer.

Note that all the states are anti-symmetric in the plane of the benzene ring. This means that they are the -orbitals of the benzene ring, i.e. they are linear combinations of the p-orbitals on each carbon atom which are perpendicular to the plane of the molecule. There are 6 such p-orbitals which means that there are 6 -orbitals.

By visually inspecting the symmetry of the different eigenstates the remaining two -orbitals are identified as state 10 and 18. These are visualized in Figure 3.

In this chapter you will analyze the transmission eigenstates at the peak positions in the transmission spectrum, i.e. at energies -3, -1.5 and 2.394 eV.

Re-open the Script Generator and delete the Eigenstate blocks.

Instead, add three Analysis/TransmissionEigenvalues blocks.

Set the energies to -3, -1.5 and 2.394 eV, into the respective TransmissionEigenvalues blocks.

Execute the script using the Job Manager and inspect the log file.

The transmission eigenvalues are obtained by diagonalizing the transmission matrix. The number of eigenvalues indicates the number of individual channels through the molecule, and the eigenvalue shows the strength of each channel. The eigenvalues are the true transmission probabilities, and thus lie in the interval [0,1]. If several channels are available at a particular energy, their sum - and hence the transmission coefficient at this energy - may however be larger than 1.

Most of the eigenvalues are very small and can be neglected. The two largest eigenvalues at each energy is listed below:

| Energy                       = -3.000000e+00 eV          |
| -------------------------------------------------------- |
| Number of transmission modes = 21                        |
+----------------------------------------------------------+
Eigenvalues (Up):
  8.108052e-01
  1.734361e-02

| Energy                       = -1.500000e+00 eV          |
| -------------------------------------------------------- |
| Number of transmission modes = 7                         |
+----------------------------------------------------------+
Eigenvalues (Up):
  7.993159e-01
  2.220198e-02

| Energy                       =  2.394000e+00 eV          |
| -------------------------------------------------------- |
| Number of transmission modes = 6                         |
+----------------------------------------------------------+
Eigenvalues (Up):
  6.740397e-01
  2.157626e-01

Note that at the first two energies there is only one dominating eigenstate, while at the last energy there are two dominating eigenvalues. (The fact that they do not sum up to the total transmission spectrum T(E) will be discussed in the next section!)

The next step is to visualize the transmission eigenstates for the first two peak energies, where there is only one significant eigenvalue.

The objects will again be saved in the file dtb_analysis.nc and can be visualized using the same approach as for the eigenstates. Below are shown isosurfaces of the two transmission eigenstates.

The transmission eigenstate at -3.0 eV has a clear resemblance with MPSH state 13, i.e. the HOMO-1 state. The transmission eigenstate at -1.5 eV is also a π-state,

In this section you will calculate the energy-dependent Local Device Density of States (LDDOS). It is related to the PDDOS calculated above, however, instead of projecting the DDOS onto orbitals, the DDOS is projected in real space.

The aim is to calculate the LDDOS in the energy range -3,3 eV with steps of 0.2 eV. This gives 31 spectra, and to set this up in the Script Generator will require that you add 30 more LocalDeviceDensityOfStates blocks. Instead of this tedious task, you will use scripting.

Use "Send To" to transfer the script to the Editor.

In the Editor, add the following lines before the LocalDeviceDensityOfStates block.

energies = numpy.linspace(-3,3,31)
for e in energies:

Indent the rest of the script by selecting the lines and pressing the TAB key. Finally, change the energy keyword into

        energy=e*eV,

The script should now look like the following

configuration = nlread('au_dtb_au.nc', object_id='gID000')[0]

# make list of energies to be used for ldos
energies = numpy.linspace(-3,3,31)

# calculate ldos for each energy in the list
for e in energies:
    local_device_density_of_states = LocalDeviceDensityOfStates(
        configuration=configuration,
        energy=e*eV,
        kpoints=MonkhorstPackGrid(3,3),
        contributions=All,
        energy_zero_parameter=AverageFermiLevel,
        infinitesimal=1e-06*eV,
        self_energy_calculator=KrylovSelfEnergy(),
        spin=Spin.Sum,
        )
    nlsave(r'lddos.nc', local_device_density_of_states)

lddos.py

The script will take around 30 min to execute on a serial machine. You may execute it using the Job Manager.

To plot the data use the script below

# import list with lddos
lddos_list = nlread('lddos.nc', LocalDeviceDensityOfStates)

#Find the z-spacing
dX, dY, dZ = lddos_list[0].volumeElement().convertTo(Ang)
dz = dZ.norm()
shape = lddos_list[0].shape()
z = dz * numpy.arange(shape[2])

# calculate average lddos along z for each spectrum
energies = []
lddos_z = []
for lddos in lddos_list:
    energies = energies + [lddos.energy().inUnitsOf(eV)]
    avg_z = numpy.apply_over_axes(numpy.mean,lddos[:,:,:],[0,1]).flatten()
    lddos_z = lddos_z + [avg_z]

# plot as contour plot
# make variables
energies = numpy.array(energies)
X, Y = numpy.meshgrid(z, energies)
Z = numpy.array(lddos_z).reshape(numpy.shape(X))


import pylab
#plot the LDDOS(E, z)
pylab.xlabel('z (Angstrom)',fontsize=12,family='sans-serif')
pylab.ylabel('Energy (eV)',fontsize=12,family='sans-serif')
contour_values = numpy.linspace(0,0.02,21)
pylab.contourf(X, Y, Z, contour_values)
pylab.colorbar()
pylab.axis([6,16,-3, 3])

pylab.title('Au-DTB-Au : LDDOS(E, z)')
pylab.savefig('lddos.png',dpi=100)


pylab.show()

lddos_ana.py

Execute the script, and you should see the following plot.

The LDDOS shows a high density in the metallic gold electrodes (white regions) and the central colored region corresponds to the vacuum gap. In the central region there are three peaks as a function of energy, positioned at -3, -1.4, 2.4 eV, in agreement with the peaks in the transmission spectrum.

If you add more LDDOS objects (i.e. more points in energy, the script above uses 31) you can enhance the energy resolution.