1. No category

A primer on Smeagol
http://www.smeagol.tcd.ie
Víctor García Suárez
Outline
1) Introduction
2) Theory
3) How to run the code
4) Simple examples
5) Some calculations
1) Introduction
Introduction to Smeagol
S pin
and
M olecular
E lectronics in an
A tomically
G enerated
O rbital
L andscape
Why Smeagol
Need for smaller electronic devices. Atomic limit
- Faster
- Cheaper
- More compact
- Similar features as today electronic elements (rectification, NDR, etc).
General features of Smeagol
- Density functional theory (DFT)
First-principles code based on localized functions: Siesta, Fireball, etc
+
- Non-equilibrium Green’s functions
Calculation of the density matrix, transmission and current under finite
biases
Smeagol Characteristics and Capabilities
Smeagol and spintronics: exploit the spin degree of freedom
- Spin polarized
- Non collinear
- Spin orbit
Calculations of both extended and isolated systems: G and k-point
calculations
Calculations of more than 100 atoms
Parallelized
2) Theory
The philosophy behind a Smeagol calculation
V
Low resistance
Large resistance
Low resistance
Direction of electronic transport (z-axis)
Left bank (reservoir)
Left electrode
Junction
Right bank (reservoir)
Extended molecule Right electrode
Leads and extended molecule
Left electrode
Molecule
Right electrode
HL [nL] – e mL
mL  m – eV/2
No m
No Fermi
distribution
HR [nR] – e mR
mR = m + eV/2
Bulk left electrode
Left reservoir
Extended molecule
Non-equilibrium
Bulk right electrode
Right reservoir
H LM
0
 H L  eV / 2



H   H ML
H M [nM ]
H MR


0
H RM
H R  eV / 2 

H M [nM ] ?
DF ( E )  f ( E )
Procedure to compute nM when
the distribution function is not Fermi

nM (r ) 
 2
ψ
(
r
 i ) f (εi )
i
Non-equilibrium
Green’s function formalism
Calculation of the leads properties
- Unit cell that is repeated along the transport direction (z)
- Use of k-points. Necessary to converge the density of states
H1
H0
H1
H0
H1
H0
H1
H1
H0
H0
H1
H0
H1
H0
H0
- At the end of the calculation:
* Hamiltonian and overlap matrices (H0, S0 and H1, S1) → Surface GF
* Density matrix
* Fermi energy
Self-energies and G matrices
The self-energies are calculated with the couplings and the surface GF
G0R is the retarded Green’s function of the leads, calculated using a semianalytic formula
0R
ˆ
ˆ  Hˆ ]
ˆR  [( E  i )Sˆ  Hˆ ]
G
(
E
)[(
E

i

)
S
L/R
M L/R L/R
L/R M
The G matrix contains information on the coupling between the extended
molecule and the leads
It is important in the calculation of the transmission
ˆG ( E)  i[ˆ R ( E)  ˆ R  ( E)]
L/R
R/L
R/L
Calculation of the surface GF
Since the GF depends on energy it is necessary to calculate k(E) from the
block vectors instead of E(k) → Solve the inverse secular equation
[ K 1e
 ik
 K 0  K1e ]k  0 ; K i  H i  ESi
ik
This involves obtaininig the inverse of the K1 matrix, which can not always be
inverted → Identify the singularities (GSVD) and get rid of them (decimation)
Once the bulk GF is constructed the surface GF is obtained by applying the
appropriate boundary conditions
The surface GF should vanish at z1 → Add a wave function to the bulk GF
z-1
z0
z1
Calculation of the extended molecule
- Includes molecule or central scattering region + part of the leads
- The first and last part of the EM must coincide with the unit cells of each of
the leads (buffer layers). This implies that:
* The same general parameters as in the bulk calculation have to be used in the
unit cell: temperature, mesh, perpendicular k-points, spin-polarization, etc
* The same particular parameters for the leads have to be used in the buffer
layers and rest of the leads in the EM: basis set, atomic coordinates, etc.
Buffer layers
Buffer layers
- The Hamiltonian, overlap and DM of the buffer layers is substituted
by those of the leads
- The buffer layers ensure that:
* The electronic structure at the beginning and end of the EM is that of
the leads
* (forced) Convergence between the electronic strucrure of the leads
near the scattering region and deep into the leads
* Absence of spurious effects in case of surfaces
TEST: infinite system
Energy, charge and
transmission
T(E)
E
Equal leads
- In case of equal leads it is possible to make the calculation of the EM
periodic to avoid the presence of surfaces. This does not mean that the
system is periodic
This does not mean that the system is periodic
No k-points are necessary
Energy mismatch between the bulk and EM calculations
In general the energy origin in a calculation of an infinite system is arbitrary
It is necessary to determine such mismatch and correct it. Otherwise the
Fermi energy can be incorrectly defined and the system can win or loose
charge
The correction is made at
certain positions of the bulk
slices
It is necessary to calculate
the bulk Hartree potential
X
X
Calculation of the DM
Siesta:
Solve eigenvalue problem. Order-N or diagonalization


*
ˆ
ˆ
Hc  ESc  ˆ   cˆn cˆn f n
n
Smeagol:
* Semi-infinite leads
* Non-equilibrium charge distribution
1


ˆ
ˆ
G ( E )  ˆ 
dE G ( E )

2i
Density matrix in equilibrium
In equilibrium it is only necessary to know the retarded Green’s function
Built with the Hamiltonian, overlap and self-energies
R
R
R
1
ˆ
ˆ
ˆ
ˆ
ˆ
G ( E )  [( E  i ) S  H   R ( E )   L ( E )]
The density matrix is obtained by integrating along an energy axis
1 
R
ˆ
ˆ
    dE[Im G ( E )] f ( E  m )
 
Complex energy contour
The lesser Green’s is calculated on an energy contour in equlibrium
Three parts: imaginary circle, imaginary line and Fermi poles





R
ˆ
ˆ
ˆ
dE G ( E )   dE G ( E )  2ik T z Im[ G ( zn )]
C
n
Fermi function poles
zn  i(2n  1)kT
Image from Atomistix
www.quantumwise.com
Brandbyge et al. Phys. Rev.
B 65, 165401 (2002)
Density matrix out of equilibrium
Out of equilibrium it the GF is not analytic inside the contour
The calculation of the DM is divided in two parts
ˆ  ˆ L/R  ˆ V
ˆ L/R
1
R
ˆ
   dE[Im G ( E )] f ( E  m L/R )
 C
1
ˆ V 
2



Rˆ
R
ˆ
ˆ
dE [G GR/LG ]( E )[ f L/R  f R/L ]( E )
Out of equilibrium voltage profile
The Hartree potential is defined up to a constant and a linear ramp (solution of
the Poisson equation)
A linear ramp related to the bias voltage is added to help the convergence out
of equilibrium
+eV/2
mL
mR
V
-eV/2
Inside the Siesta version of Smeagol
Initial guess
Calculate effective potential
Smeagol
Calculate the DM using NEGF
ˆ 
1
dE Gˆ  ( E )

2i
Smeagol
Smeagol
No
Compute electron density
Yes
Self-consistent?
Output quantities
Transmission and
current
Calculation of the transmission and current
The transmission is calculated at the end of the self-consistent cycle
It is possible to simplify it due to the small size of the G matrices
Rˆ ˆR
ˆ
ˆ
T ( E )  Tr [GLG GR G ]( E ) 
Rˆ ˆR
ˆ
ˆ
 Tr [G G G G ]
L
LR
R
RL
The current is calculated by integrating the transmission
e 
I   dE T ( E )[ f ( E  m L )  f ( E  m R )]( E )
h 
Single level coupled to wide band leads
Wide band leads have a constant density of states at the Fermi level
e1
After coupling the level the onsite energy
(e0) is renormalized by the real part of the
self-energy (e1)
The imaginary part of the self energy give
the inverse of the lifetime (width of the
Breit-Wigner resonance)
ˆ RL/R    i G
2
2
e1  e 0  
2
1
G
G R (E) 
T (E) 
E  e 1  iG
( E  e1 ) 2  G 2
3) How to run the code
Calculation of the leads
First run before calculating the transport properties of the extended molecule
- Include two variables in the input file:
BulkTransport T . It specifies wehter or nor the bulk parameters are written
BulkLead
LR . Left (L), right (R) or both (LR) leads
* At the end of the calculation three or four files are generated:
- bulklft.DAT and bulkrgt.DAT: contain the label of the system and basic
information such as the Fermi energy, temperature, etc.
- SystemLabel.HST: contains the Hamiltonian and overlap matrices
- SystemLabel.DM: contains the density matrix
Example of bulk file
SystemName
Au
SystemLabel
Au
NumberOfAtoms
2
NumberOfSpecies
1
%block ChemicalSpeciesLabel
1 79 Au
%endblock ChemicalSpeciesLabel
%block PAO.Basis
Au 1
n=6 0 1
5.0
%endblock PAO.Basis
%block Ps.lmax
Au 1
%endblock Ps.lmax
LatticeConstant
1.00 Ang
%block LatticeVectors
10.00 0.00 0.00
0.00 10.00 0.00
0.00 0.00 4.08
%endblock LatticeVectors
AtomicCoordinatesFormat Ang
%block AtomicCoordinatesAndAtomicSpecies
0.00 0.00 0.00 1 Au 1
0.00 0.00 2.04 1 Au 2
%endblock AtomicCoordinatesAndAtomicSpecies
%block kgrid_Monkhorst_Pack
1 0 0 0.0
0 1 0 0.0
0 0 100 0.0
%endblock kgrid_Monkhorst_Pack
xc.functional
GGA
xc.authors
PBE
MeshCutoff
200. Ry
MaxSCFIterations
10000
DM.MixingWeight
0.1
DM.NumberPulay
8
DM.Tolerance
1.d-4
SolutionMethod
diagon
ElectronicTemperature
150 K
SaveElectrostaticPotential
T
BandLinesScale pi/a
%block BandLines
1 0.00 0.00 0.00
200 0.00 0.00 0.50
%endblock Bandlines
BulkTransport
T
BulkLead
LR
DM.UseSaveDM
T
Calculation of the extended molecule. General
Run the extended molecule file in the same directory that contains the bulk files
- Variable to define the transport calculation:
EMTransport
T . Performs the transport calculation
- Variables related to the energy contour:
NEnergReal
500 . Number of points along the real axis (out of equilibrium)
NEnergImCircle 90 . Number of points in the imaginary circle
NEnergImLine
30 . Number of points in the imaginary line
Npoles
10 . Number of poles of the Fermi function
Delta
1.0d-4 . Small imaginary part of the Green’s Function
EnergyLowestBound -10.0 Ry . Energy of the lowest bound of the EC
Nslices
1 . Number of slices that are substituted by the bulk H and S
Calculation of the extended molecule. Out of equilibrium
- Variables related to the bias voltage (out of equilibrium):
VInitial
0.0 . Initial value of the bias voltage
VFinal
2.0 . Final value of the bias voltage
NIVPoints
AtomLeftVcte
AtomRightVcte
10 . Number of points where the bias is going to be applied
9 . Left position where the bias ramp starts
36 . Right position where the bias ramp ends
At each bias point the electronic structure is converged. Optionally, it is also
possible to relax the atomic coordinates
When the electronic structure converges the transmission and current at that
bias point are calculated
Calculation of the extended molecule. Transmission and VH
- Variables related to the calculation of the transmission
NTransmPoints 800 . Number of energy points where the transmission is
going to be calculated
InitTransmRange
-8.0 eV . Initial value of the transmission range
FinalTransmRange 2.0 eV . Final value of the transmission range
- Variables related to the energy mismatch between bulk and EM
HartreeLeadsLeft
HartreeLeadsRight
2.040 Ang . Left position where the correction is applied
10.200 Ang . Right position where the correction is applied
HartreeLeadsBottom -1.711 eV . Value of the Hartree potential at a certain point in
the leads
Example of extended molecule file
SystemName
Au.em
SystemLabel
Au.em
NumberOfAtoms
6
NumberOfSpecies
1
%block ChemicalSpeciesLabel
1 79 Au
%endblock ChemicalSpeciesLabel
%block PAO.Basis
Au 1
n=6 0 1
5.0
%endblock PAO.Basis
%block Ps.lmax
Au 1
%endblock Ps.lmax
LatticeConstant
1.00 Ang
%block LatticeVectors
10.00 0.00 0.00
0.00 10.00 0.00
0.00 0.00 12.24
%endblock LatticeVectors
AtomicCoordinatesFormat Ang
%block AtomicCoordinatesAndAtomicSpecies
0.00 0.00 0.00 1 Au 1
0.00 0.00 2.04 1 Au 2
0.00 0.00 4.08 1 Au 3
0.00 0.00 6.12 1 Au 4
0.00 0.00 8.16 1 Au 5
0.00 0.00 10.20 1 Au 6
%endblock AtomicCoordinatesAndAtomicSpecies
...
EMTransport
T
NEnergReal
500
NEnergImCircle
60
NEnergImLine
30
NPoles
10
Delta
2.d-4
EnergLowestBound -8.d0 Ry
NSlices
1
VInitial
0.d0 eV
VFinal
2.d0 eV
NIVPoints
10
AtomLeftVCte
2
AtomRightVCte
7
TrCoefficients
T
NTransmPoints
800
InitTransmRange -14.0 eV
FinalTransmRange
8.0 eV
PeriodicTransp
T
UseLeadsGF
F
HartreeLeadsLeft
2.040 Ang
HartreeLeadsRight 10.200 Ang
HartreeLeadsBottom -1.711 eV
#%block SaveBiasSteps
#012
#%endblock SaveBiasSteps
DM.UseSaveDM
T
4) Simple examples
Infinite atomic chain. Equilibrium
Band structure and transmission
Two atoms in the unit cell → two
crossing bands in the Brillouin zone
One single channel at every energy
Perfect squared transmission
Infinite atomic chain. Out of equilibrium
Out of equiilibrium transmission
Starts to disappear at the edges, where
the bands of both leads mismatch
Current
Ohmic regime
Diatomic molecule. Equilibrium I
Effect of separating the atoms from the leads
Changing the coupling
configuration
Diatomic molecule. Equilibrium II
Effect of changing the intramolecular distance
Decrease the distance
Increase the distance
Diatomic molecule. Out of equilibrium
Bonding
Antibonding
Different behaviour of the bonding and antibonding orbitals
Bias-dependent transmission
Negative differential resistance (NDR)
4) Some calculations
Magnetoresistance effects in nickel chains
Properties of this system:
- Highest magnetic moments in the
middle of the chain
- The spins invert close to the leads
Spin-polarized
Non-collinear
- Abrupt change (collinear) of the
magnetization
Symmetric, parallel
Symmetric, antiparallel
Aymmetric, antiparallel
Metallocenes inside carbon nanotubes
Chains of metallocenes inside CNTs
Reduction of the magnetoresistance due to
charge transfer from the metallocene
CoCp
CoCp +
CNT
CNT
Magnetoresistance effect
Even-odd effect in monoatomic chains
Different conductance depending on the
number of atoms in the chain
The oscillations depend on the type of contact
configuration
Gold
Sodium
Conductance of H2 molecules couples to Pt or Pd leads
Two possible configurations: parallel or perpendicular to the transport direction
In case of Pd the H2 molecule can go inside the bulk and contacts
Continuous
line: Pt
Dashed line:
Pd
Oscillations in Pt chains
Chains with a number of atoms between 1 and 5
Structural oscillations due to changes in the levels at the Fermi level as the
chain is stretched from a zigzag to a linear configuration
2 atoms
3 atoms
I-V calculations of polyynes between gold leads
Two types of molecules with different
coupling atoms
S
Two types of NDR due to a different
evolution of the resonances with bias
N
(a)
(b)
Porphyrins between gold leads
Very large calculations with more than 500 atoms
Evolution of the conductance with the number of porphyrin units and the angle
between them
Fin