In [ ]:
from __future__ import print_function
import sisl
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
$ \newcommand\HH{\mathbf{H}} \newcommand\SO{\mathbf{S}} \newcommand\Scat{\boldsymbol\Gamma} \newcommand\ket[1]{|#1\rangle} \newcommand\bra[1]{\langle#1|} \newcommand\set[1]{\{#1\}} \newcommand\eig{\epsilon} $ In this example you will perform transport calculations by projecting scattering states onto molecular orbitals.
We will select Carbon chains as our electrodes and the C$_{60}$ fullerene as our molecule (it has nice symmetries). The molecular projection method are created from a subset molecule space ($\set M$): \begin{equation} \HH_{\set M} \ket{i} = \eig_i^{\set M}\SO_{\set M} \ket{i} \end{equation} with the projectors orthogonalized through the Löwdin transformation \begin{equation} \ket{i'}=\SO^{1/2}\ket i \end{equation} Then the projection of the scattering states will read \begin{equation} \tilde\Scat = \sum\ket{i'}\bra{i'}\Scat \sum \ket{j'}\bra{j'} \end{equation} which completes the basis transformation.
These projectors can be performed in TBtrans by defining the $\set M$ region and defining which scattering matrices should be projected onto.
From the above few equations it should be obvious that to create such a projection it is necessary to define a device region where the scattering matrices are living only on the projected region (i.e. the molecule).
Instead of manually setting up the C$_{60}$ molecule we find the coordinates on this web-page (http://www.nanotube.msu.edu/fullerene/fullerene-isomers.html). The coordinates are stored in the C60.xyz file in this directory. Note that when reading this geometry it does not know the orbital distance so we have to calculate it.
In [ ]:
C60 = sisl.Geometry.read('C60.xyz')
# Calculate the nearest neighbour distance
dist = C60.distance(R=5)
C60.atom.atom[0] = sisl.Atom(6, R=dist[0] + 0.01)
print(C60)
Now we have the molecule, all we need are some electrodes and the connection from the electrodes to the molecule. This electrode is a 2x2 square lattice with transport along the $x$-direction.
In [ ]:
elec = sisl.Geometry([0] * 3, sisl.Atom(6, R=1.), [1, 1, 10]).tile(2, 0).tile(2, 1)
elec.set_nsc(a=3, b=1)
H_elec = sisl.Hamiltonian(elec)
H_elec.construct(([0.1, 1.1], [0., -1]))
H_elec.write('ELEC.nc')
Create the final device by making the electrode have 2 screening layers, then the C$_{60}$ and finally the right-electrode (equivalently setup to the left part).
In [ ]:
elec_x = elec.tile(3, 0)
# Translate to origo
C60 = C60.translate(-C60.center(what='xyz'))
C60 = C60.translate([-C60.xyz[:, 0].min(), 0, 0])
# Do trickery to make sure the coordinates are consecutive along x
C60.set_sc([C60.xyz[:, 0].max() + 1., 1., 1.])
device = elec_x.append(C60, 0).append(elec_x, 0)
The full device is now created and we simply need to create the electronic structure.
In [ ]:
H = sisl.Hamiltonian(device)
H.construct(([0.1, 1.1], [0., -1]))
# Correct the C_60 couplings to something different
idx_C60 = np.arange(len(elec_x), len(elec_x) + len(C60), dtype=np.int32)
for ia in idx_C60:
idx = device.close(ia, R=[0.1, C60.maxR()])[1]
# On-site is already 0, so don't bother doing anything there
# Split idx into C60 couplings and chain couplings
for i in idx:
if i in idx_C60:
H[ia, i] = -1.5
else:
# Coupling to chain
# Since we are only looping atoms in C60
# we have to also set the coupling into C60
# (to assert Hermiticity)
H[ia, i] = 0.1
H[i, ia] = 0.1
H.write('DEVICE.nc')
Calculating projceted transmissions is verbose because multiple things is required to be defined. The following flags are required to be specified, careful read about each of them:
C60)